基于集中质量、仿真和灰色关联分析法对系泊系统的研究

基于集中质量、仿真和灰色关联分析法对系泊系统的研究 摘要

系泊系统是近浅海观测网的重要组成部分，如何设计最优的系泊系统从而促 进信号的传播具有重要的意义。

本文针对最优系泊系统的设计进行研究，主要解决了如何确定锚链的型号、

长度和重物球的质量使得浮标的吃水深度和游动区域及钢桶的倾斜角度尽可能 小的问题，同时考虑在不同风速，不同海水速度下的锚链形状和浮标的游动区域 的问题，给出了最优系泊系统设计的方案。

针对问题一我们采用集中质量法将系统中的各个物体视为一个质点，对各个 物体建立静力平衡方程，在水深 18m 时给定浮标在海水中所受浮力，从而根据我 们建立的平衡方程求出各物体的倾斜角度，再根据几何关系求出海域的模拟深度，

通过不断修正浮标的浮力使得海域的模拟深度等于 18m，最终求得风速分别为 12m/s 和 24m/s 时浮标的吃水深度 h 为 0.7397m 和 0.74883m，同时给出不同风

速下钢桶和各节钢管的倾斜角度及浮标的最远位置（浮标的游动区域视为一个圆 面）。考虑到锚链由 210 节链环构成，通过对每节链环进行受力分析确定了每节 链环的位置，从而给出了链环的形状图像。在求解过程中由于拉力具有不确定性，

我们通过两次角度代换使得程序可以顺利地运行。

针对问题二我们沿用了问题一的算法，求得风速为 36m/s 时钢桶和各节钢管 的倾斜角度，锚链形状和浮标的游动区域。其中钢管的倾斜角度   7.998

 5

， 锚链末端与锚的链接处的切线方向与河床的夹角   17.829

 16

，不符合设备 的工作要求。通过绘制不同质量的重物球分别和   , , h

的散点图，得出了重物球 的质量与 ,   成负相关，与 h 成正相关的结论，同时绘制了不同质量的重物球下

系泊系统的大致图像。然后沿用问题一的算法，在   5 ,   16 , h  2 的限制条 件下，逐渐增加重物桶的质量，求得重物球的质量范围为  ^{1800 kg ,5200 kg}  _。

针对问题三的第一个子问题我们首先考虑最恶劣的环境，即布放点的海水速 度达到最大值 1.5m/s、风速达到最大值 36m/s。接着利用穷举法给出了锚链型号、

长度和重物球质量不同组合的 1000 种系泊系统设计方案，在问题 2 的限制条件 下遴选可行的方案，在不同的水深下建立灰色关联分析综合评价模型，计算不同 水深的关联度值，将同一方案的关联值相加，最终将关联度值最高的设计方案确 定为最优的系泊系统的设计，得到了最优的系泊系统设计方案为锚链的型号为V ， 长度为 20m，重物球的质量为 3900kg。针对问题三的第二个子问题，我们在最优 的设计方案下，改变风速和海速，得到了不同情况下钢桶、钢管的倾斜角度、锚 链形状、浮标的吃水深度和游动区域。

关键词：系泊系统；平衡受力分析；集中质量法；仿真；灰色关联分析法

1、 问题重述

本题着重讨论了系泊系统的设计，以相关的物理知识为基础，利用几何求解 和最优规划等数学工具，解决以下几个问题：

问题一：在电焊锚链的型号和长度、重物球质量、水深、海水密度已知的情 况下，不考虑海水的流动，分别计算海面风速为 12m/s、24m/s 时钢桶和各节钢 管的倾斜角度、锚链形状、浮标的吃水深度和游动区域。

问题二：问题一的假设不变，风速为 36m/s 时可能会出现钢桶的倾斜角超过 5 度，锚链末端与锚的链接处的切线方向与河床的夹角超过 16 度的现象，如何 调节重物球的质量，使得系泊系统符合设计的角度要求。

问题三：将问题三分解为两个子问题

1.考虑潮汐等因素的影响，如何设计系泊系统使得浮标的吃水深度和游动区 域及钢桶的倾斜角度尽可能小。

2：布放海域的风力、水流力和水深情况在一定的范围内变动，分析不同情 况下钢桶、钢管的倾斜角度、锚链形状、浮标的吃水深度和游动区域。

2、问题分析

问题一的分析：

考虑到钢桶和各节钢管的倾斜角度、锚链形状、浮标的吃水深度和游动区域 与整个系泊系统的各部分受力情况息息相关，于是对浮标、各节钢管、钢桶和各 个锚链一一进行受力分析，在浮标、钢管、钢桶、锚链和锚平衡状态下，得出它 们的静力平衡方程，从而计算出各拉力和各个倾斜角度，由于钢管、钢桶、锚链 长度已知，结合角度关系，建立以锚链的末端与锚的链接点为坐标原点的直角坐 标系，通过 matlab 编程绘制出整个系泊系统的图像，最后通过几何关系判断锚 是否会着底。在满足水深为 18 m 的条件下，求得海面风速分别为 12 m s 和 24 /

/

m s 时钢桶和各节钢管的倾斜角度、锚链形状、浮标吃水深度和游动区域。其中 游动区域是以游动区域半径衡量的，游动区域半径长度指的是锚向海平面的竖直 投影点到浮标平衡位置的距离，游动区域是一个圆面。

问题二的分析：

问题一的条件不变，沿用问题一的算法可以计算出海面风速为 36 m s 时钢 / 桶和各节钢管的倾斜角度、锚链形状、浮标的吃水深度和游动区域；

考虑到在重物球为 1200 kg 时，钢桶的倾斜角超过了 5 度，锚链末端与锚的 链接处的切线方向与河床的夹角超过了 16 度，不改变其他条件分析可知重物球 的质量与钢桶的倾斜角  和锚链末端与锚的链接处的切线方向与河床的夹角  成负相关，即重物球质量越大，  和  越小，故在   5 ,

  16 ,

h

 2 条件下，

可以通过不断调节重物球的质量，找到重物球的最小质量和最大质量。

问题三的分析：

对问题三第一个子问题的分析：

系泊系统的设计问题就是要确定锚链的型号、长度和重物球的质量，使得浮 标的吃水深度和游动区域及钢桶的倾斜角度尽可能小。

沿用问题一的算法，增加一个与风力方向相同的水流力，近海水水流力可以

通过近似公式 F

=374  S

 V

计算。考虑最恶劣的环境，布放点的海水速度达 到最大值 1.5m/s、风速达到最大值 36m/s。根据穷举法原理，给出锚链型号、长 度和重物球质量不同组合下的 1000 种系泊系统的设计方案，水深在 16 米至 20 米的区间内变化，我们可以从中遴选出钢桶倾斜角度小于 5 度，锚链末端与锚的 链接处的切线方向与河床的夹角小于 16 度的设计方案，在不同的水深下建立灰 色关联分析综合评价模型，计算不同水深的关联度值，将同一方案的关联值相加，

最终将关联度值最高的设计方案确定为最优的系泊系统的设计。

对问题三第二个子问题的分析：

在子问题一所确定的最优系泊系统设计的条件下，可以通过不断改变海水速 度和风速，分析不同情况下钢桶、钢管的倾斜角度、锚链形状、浮标的吃水深度 和游动区域。

3、模型假设

3.1 假设钢管与钢管链接点的质量忽略不计。

3.2 假设重物球和锚链在海水中所受浮力忽略不计。

3.3 假设海面风速和海水流速为单向流，且方向一致沿 _X 轴正方向。

3.4 假设近浅海观测网的传输节点位于同一平面。

3.5 假设在海水流动情况下，仅考虑水流力对浮标、钢管和钢桶的影响，对 锚链和重物球的影响忽略不计。

3.6 假设在海水流动情况下，锚受到的水流力不影响锚链末端与锚的链接处 的切线方向与海床的夹角。

3.7 假设链环在拉力作用下形状不发生改变。

4、符号说明

4.1 符号说明

浮标的符号说明：

符号 符号说明 符号 符号说明

m

浮标质量 h

浮标的吃水深度

g 重力加速度 F

浮标所受海水浮力

G

浮标重力 S

浮标在风向法平面 的投影面积



海水密度 V

风速

V

浮标排开海水的体积 W

浮标受到的海风载 荷力

s 圆柱体底面积 T

第一节钢管对浮标

的拉力

钢管的符号说明：

符号 符号说明 符号 符号说明



即第 i 节钢管倾斜角度 T

第 i  1 节钢管对第 i 节钢管的拉力 m

每节钢管质量 T

第 i 节钢管对第 i  1

节钢管的拉力 G

每节钢管重力 T

浮标对第一节钢管

的拉力 d

每节钢管直径 T

钢桶对第四节钢管

的拉力 h

每节钢管长度 T

第四节钢管对钢桶

的拉力 F

每节钢管所受海水浮力 V

_{钢管的体积}

钢桶的符号说明：

符号 符号说明 符号 符号说明

m

钢桶质量 F

钢桶所受海水浮力 G

钢桶重力 h

钢桶的长度 m

重物球质量 d

钢桶的外径

T

重物球拉力  钢桶的倾斜角度

链环（锚链）的符号说明：

符号 符号说明 符号 符号说明

m

每节链环的质量 L

锚链的长度

G

每节链环的重力 p

^{锚链为 j 型号时的单}

位长度质量

'

T

第一节链环对钢桶拉力 T

第 j+1 节链环对第 j 节链环的拉力

'

T

钢桶对第一节链环的拉力 T

第 j 节链环对第 j+1 节链环的拉力



第 i 节链环的倾斜角度  ^{锚的链接处的切线方}

向与海床的夹角

问题三的符号说明

符号 符号说明 符号 符号说明

F 近海水水流力  (k) 比较数列 x 在第 k

个指标的关联系数

S

物体在水流速度法平面的投

影面积  _分辨系数

V

水流速度 ^{min min}

^x

^{(t)  x}

^(t) 两级最小差

 _权重 max max

(t)

(t)

s t

x  x 两级最大差

x

参考数列 r

第 i 个设计方案的

灰色关联度值

5、模型建立与求解

5.1 问题一：在海水静止时不同风速下系泊系统的数学模型 求解过程的流程图如下：

对浮标进行设计分析 假设浮标的浮力 在浮标受力平衡下求出对

一号钢管的拉力与其夹角

依次对各钢管、钢桶、锚 链中的各个链环进行受力 分析，迭代求出各物体的

拉力及其倾斜角

将物体高度乘以倾斜角的 余弦值得出系泊系统竖直

方向的高度 判断锚是否着地（即计算

的高度与题目所给的海底 深度之间是否相等）

误差范围外

（修正浮标的浮力）

误差范围内 记录所需的倾斜夹角、吃

水深度、游动区域并用程 序仿真锚链形状、

5.1.1 系泊系统求解分析和模型建立

考虑静止情况下的集中质量受力分析法，注意到海水静止，故不需考虑海水 对钢管、钢桶和锚链的水流阻力。现分别对系泊系统中的物体进行受力分析：

浮标受力分析

设风力方向沿 _X 轴正方向，为单向平面流。在海水静止时，浮标受到的力有

重力 G

，浮力 F

,海风荷载力 W

和钢管拉力 T

，如图 1 所示

图 1：浮标静力平衡状态的受力示意图

其中 x 轴正方向为锚向海平面的竖直投影点指向浮标，Y 轴正方向为竖直向 上，后面图形中的方向取向与图 1 相同。

根据浮标质量，计算浮标重力：

0 0

G  m g （1-1）

浮标所受浮力由其排开海水的体积和海水密度决定，浮标为圆柱体，所以其 浮力计算为：

0

=

F  

gV



gsh （1-2）

式(1-2)中



表示海水密度；V

表示浮标排开海水的体积； ^s 表示圆柱体底面积； h

表示浮标的吃水深度

由说明可知，近海风荷载可通过近似公式 W

=0.625 SV 计算：

2 0

=0.625

W S V 、 S

   2 (2 h

) （1-3）

式（1-3）中，

S

表示浮标在风向法平面的投影面积； V

表示风速

浮标处于平衡状态，将各个分力正交分解，静力平衡方程由（1-4）式给出：

0 0 0 1

0 0 1

cos sin

F G T

W T





 

 竖直方向上：

水平方向上： （1-4）

式（1-4）中，



表示第一节钢管的倾斜角度

求解（1-4）静力平衡方程，可以得出浮标下端钢管拉力 T

和第一节钢管的

倾斜角度 

：

2 2

0 0 0 0

0 1

0 0

( )

arctan

T F G W

W F G



   

 

  



（1-5）

在上述计算中，由于浮标在海水中所受浮力是未知的，从而导致浮标的吃水 深度无法计算得出，需要给定浮力初始值 F

，从而得到 T

和 

。

钢管受力分析

将钢管看作一个质点，对每节钢管一一受力分析，如图 2 所示，由于每节钢 管规格一致，即质量、长度和直径相同，所有每节钢管重力和所受浮力相等。每 节钢管受到的力有钢管自身重力 G

、浮力 F

、上一节钢管拉力 T

和下一节钢 管拉力 T

，( i  1, 2, 3, 4) ，其中 T

表示浮标对第一节钢管的拉力， T

表示钢桶 对第四节钢管的拉力。

图 2：钢管静力平衡状态的受力示意图 根据牛顿第三定律，力的作用相互性可知，

0 1,0

, 1 1,

i i i i

T T T

T

 

 

 （2-1）

根据每节钢管质量，计算每节钢管质量：

gg gg

G  m g （2-2）

已知每节钢管直径和长度的条件下，根据浮力公式计算出每节钢管的海下 浮力：

( )

4

gg

gg gg gg

F  

gv  

g pid h （2-3）

每节钢管处于平衡状态，将各个分力正交分解，静力平衡方程由（2-4）式 给出：

竖直方向： G

 T

cos 

 F

 T

cos 

水平方向上： T

sin 

 T

sin 

（2-4）

式 （2-4）中，



表示第 i 节钢管的倾斜角度

在已知 T

、 

、 G

、 F

的前提条件下，由（2-4）静力平衡方程可以求得 第二节钢管对第一节钢管的拉力 T

及第二节钢管的倾斜角度 

，依次迭代下去，

可以求得各节钢管对上一节钢管的拉力 T

和各节钢管的倾斜角度 

：

2 2

, 1 , 1 , 1

, 1 1

, 1

( cos ) ( sin )

arctan sin

cos

i i gg i i i gg i i i

i i i

i

gg i i i gg

T F T G T

T

F T G

 

 



  







      

  

    



（2-5）

在 matlab 实现迭代的程序中，由于下一个的拉力具有不确定性，我们规定



为从 Y 正向顺时针旋转到下一个拉力方向的角度，如图 3：

图 3：matlab 迭代中角度的规定取向

由于 matlab 的反正切函数根据（2-5）求出来的角度范围为 , 2 2

    

 

  ，与我

们规定的角度取向不同，故我们做以下四种情形做四种代换：

图 4.1： T

 0, T

 0

图 4.2： T

 0, T

 0

图 4.3： T

 0, T

 0 图 4.4： T

 0, T

 0

图示中的 T

， T

分别代表钢管受到上一个钢管的拉力，自身的重力，自身 所受到的浮力的合力水平方向的分力和竖直方向的分力， T 代表钢管受到的下 一个钢管的拉力。由（2-5）可知 T

 T

 sin , 

T

 F

 T

 cos 

 G

.

记第 i 个钢管分别受第 i -1，第 i +1 个钢管拉力真实角度值为 

， 

， matlab 求解出的值为 

， 

。针对上图四种情形做的代换分别为

* '

* '

* '

* '

2 (4.1)

(4.2) (4.3) (4.4)

i i

i i

i i

i i

  

  

  

 

  

  

  

 



也有可能出现下列情形：

（i） T

 0

 T

 0  T

 0 此时 

   

此时 

 

（ii） T

 0

 T

 0  T

 0

此时 ³

i

2

   此时

i

2

   （  

,

^{的关系类似可}

得）

作出两个钢桶直接的受力分析图，如图 5 所示：

图 5：两个钢管之间的受力分析图示 为了利于 matlab 的迭代，我们再继续进行下列两种代换：

（i）

0  

  （ii）   

 2  此时 

 

  此时 

 

  钢桶受力分析

钢桶竖直时，水声通讯设备的工作效果最佳。若钢桶倾斜，则影响设备的 工作效果，为了控制钢桶的倾斜角度，钢桶与电焊锚链链接处悬挂了重物球，

考虑到重物球体积无法确定且重物球在海水中所受浮力相对其自身重力较小，

所以忽略重物球的浮力。如图所示，钢管受到的力有自身重力 G

、重物球拉力

T

、浮力 F

、第四节钢管拉力 T

和锚链拉力 T

。

图 6：钢桶静力平衡状态的受力示意图 根据钢桶质量，计算钢桶重力：

gt gt

G  m g （3-1）

由于忽略重物球所受浮力，则重物球对钢桶拉力即重物球自身重力：

zwq zwq zwq

T  G  m g （3-2）

已知钢桶外径和长度的条件下，根据浮力公式计算出钢桶的海下浮力：

2

( ) 4

gt

gt gt gt

F  

gv   

g d h （3-3）

钢桶处于平衡状态，将各个分力正交分解，静力平衡方程由（3-4）式给 出：

5,4 4

5,4 4 4,5 4

cos cos

sin sin sin

gt zwq ml gt

ml

G T T F T

T T T

 

  

   

 

竖直方向上：

水平方向上： （3-4）

在已知 G

、 T

、 F

和 T

的基础上，由（3-4）平衡方程可以得出锚链 拉力 T

和钢桶的倾斜角度  。

2 2

5,4 4 4,5 4

4,5 4

5,4 4

( cos ) +(T sin ) T sin

arctan

cos

ml gt gt zwq

gt gt zwq

T F T G T

F T G T

 

 



    

 

    



（3-5）

链环受力分析

锚链的受力分析如图 7 所示，锚链受到的力有锚链重力 G

、锚的拉力 T

和钢桶拉力 T

。由于锚链的长度过长（22.05m），这是一个非共点力的平衡问

题，不适合再采用共点力平衡的分析方法，且考虑到锚链的各个链环有着不同

的受力状态，故对锚链的 210 个 2 型链环依次进行受力分析，将每个链环视为

一个质量集中点。考虑到锚链直径题目中未给出，且搜集数据发现锚链直径非

常小，导致锚链在海水受到的浮力太小故忽略锚链在海水中所受浮力，采用和 上面钢管类似的分析方法，受力示意图如图 8：

图 7：锚链静力平衡状态的受力示意图 图 8：链环静力平衡状态的受力示意图 根据链环长度及单位长度质量，计算链环重力：

lh lh lh j

G  m g  L   p g _（4-1）

式 （4-1）中 L

表示链环的长度； p

表示锚链为 j 型号时的单位长度质量 根据牛顿第三定律，力的作用相互性可知，

' '

1,0 0,1

' '

, 1 1,

i i i i

T T T

T

 

 

  （4-2）

式 （4-2）中，

' , 1

T

表示第 j+1 节链环对第 j 节链环的拉力； T

表示第 j 节链环对第 j+1 节链环的拉力

链环处于平衡状态，将各个分力正交分解，静力平衡方程由（4-3）式给 出：

T

' '

, 1 1 , 1

' '

, 1 1 , 1

cos cos sin sin

lh i i i i i i

i i i i i i

G T T

T

 

 

  

  

   

   _（4-3）

式（4-3）中，



表示第 i 节链环的倾斜角度

在已知 T

、 

、 G

的前提条件下，由（4-3）静力平衡方程可以求得第二

节链环对第一节链环的拉力 T

及第二节链环的倾斜角度 

，依次迭代下去，可

以求得各节链环对上一节钢管的拉力 T

和各节链环的倾斜角度 

：

' ' 2 ' 2

, 1 , 1 , 1

' , 1

1 '

, 1

( cos ) ( sin )

arctan sin

cos

i i i i i lh i i i

i i i

i

i i i lh

T T G T

T

T G

 

 



  







     

  

 

 



（4-4）

5.1.2 海水静止情况下的静力学模型的构建与求解

由上文受力分析可知，在给定浮标浮力的初始值条件下，可以通过 matlab 编程计算在海水静止时海面不同风速下的钢桶和各节钢管的倾斜角度、锚链形 状、浮标的吃水深度和游动区域。因此可以通过倾斜角度的余弦值与钢桶和各 节钢管的长度的乘积计算出水深，由于水深为固定值（已知），因此可以通过 不断地修正浮标浮力的值使得计算出的水深为给定值，也即锚恰好着底。

Matlab 程序见附件 1，求得的结果如表 1 所示：

表 1：风速分别为 12m/s,24m/s，36m/s 时的求解结果

12m/s 的风速 24m/s 的风速 36m/s 的风速 第一节钢管倾斜角度 0.959656 3.726518 7.826934 第二节钢管倾斜角度 0.965234 3.747668 7.868896 第三节钢管倾斜角度 0.970877 3.769059 7.911307 第四节钢管倾斜角度 0.976587 3.790696 7.954174 钢管的倾斜角度 0.982364 3.812581 7.997506 吃水深度 0.739700 0.748830 0.769920 游动区域半径 11.615149 17.420198 18.712474 锚点与海床的夹角 0.000000 0.000000 17.828720 实际算出的水深 18.001377 18.003682 18.001244 浮标的移动区域视为以锚向海平面的投影点为中心的圆面，圆的边界位置 为浮标的平衡位置，其半径即为表 1 给出的游动区域半径。

风速分别为 12m/s,24m/s 时锚链的形状如下图 9 所示：

图 9：风速分别为 12m/s 和 24m/s 时系泊系统形状的预分析

从图中可以看出锚链的形状出现了反常的“上翘”情况，这是因为实际情 况下锚链到达海底会受到海底的支持力，而程序中没有考虑这种情况，因此在 静力平衡状态下，程序计算出的链环会给予一个沿 X 轴负方向斜向上的拉力，

故在 12m/s,24m/s 风速下 X 锚链的尾端有一段平躺于海底，平躺的链环个数分

别为 85 和 5，β 均等于 0°。调整程序代码，最终得到实际锚链的形状如图 10 所示。

图 10：风速分别为 12m/s 和 24m/s 时系泊系统形状 5.2 问题二：海水静止时不同重物球质量下系泊系统的数学模型 5.2.1 风速为 36m/s 时的模型求解

和问题一的情形一样，海面的风速变为 36m/s,带入问题一的程序中可求得 钢桶和各节钢管的倾斜角度，浮标的吃水深度和游动区域，结果见表 1。锚链 的形状如下图 11 所示：

图 11：在风速在风速为 36/s 时系泊系统的形状 5.2.2 重物球的质量区间的求解

值得注意的是，在上节中我们求得  =7.998°>5°，  =17.829°>16°，

表明在此情况下重物球的质量已经不能满足题目的要求。考虑到重物球的质量

与钢桶的倾斜角  和锚链在锚点与海床的夹角  成负相关，即重物球质量越

大，  和  值越小，我们取重物球的质量从 1200kg 开始逐渐增加它的质量到

10000kg（步长为 100 kg ），在   5 ,   16 , h  2 的条件下，得到重物球的质

量范围 1800kg-5200kg，不同质量的重物球下系泊系统的形状见图 12，matlab 求解程序见附件 2。

图 12：不同质量的重物球下系泊系统的形状 5.3 问题三: 不同情况下评判与分析最优系泊系统的数学模型 5.3.1 第一个子问题的求解：

采用穷举法，设计出 1000 种系泊系统的设计方案。设计方案包括确定锚链 的型号、长度和重物球的质量，其中锚链的型号分为 5 种；根据水深范围 16 米 至 20 米以及钢管和钢桶总长度为 5 米，确定锚链的长度范围为 16 米至 25 米，

设定 10 种锚链长度（步长为 1m）；由问题二可知调节后的重物球质量区间，

又问题三中考虑水流力的影响，故确定重物球的质量范围为 2000 kg 至 3900 kg ,设定 20 种重物球的质量（步长为 100kg）。建立灰色关联分析综合评价模 型，通过 MATLAB 编程，根据关联度值大小排名确定最终的系泊最优设计，即确 定锚链的型号、长度和重物球的质量。

灰色关联分析综合评价法模型的建立：

一、灰色关联分析法的实现步骤：

step 1、确定评价指标和评价标准（参考数列）。

系泊系统的设计目标是使得浮标的吃水深度和游动区域及钢桶的倾斜角度 尽可能小，因此评价指标为吃水深度、游动区域面积和钢桶的倾斜角度。

step 2、确定各指标值对应的权重 ^ ：

假设吃水深度、游动区域半径和钢桶的倾斜角度对系泊系统的设计影响重 要性相同，因此确定 =( , , ) ^{1 1 1}

3 3 3

 。

step 3、计算灰色关联系数：

0 0

0 0

min min ( ) ( ) max max ( ) ( ) ( ) ( ) ( ) max max ( ) ( )

s s

s t s t

i

i s

s t

x t x t x t x t

k x t x t x t x t

 



  

    （5-1）

为比较数列 x

对参考数列 x

在第 k 个指标上的关联系数，其中 ^{ }   ^{0,1 为}

分辨系数，其中称 min min

( )

( ) max max

( )

( )

s t

x t  x t 、

x t  x t 分别为两级最小差 及两级最大差。分辨系数  越大，分辨率越大；  越小，分辨率越小。

step 4、计算灰色加权关联度：

灰色加权关联度计算公式为

1

( )

n

i i i

k

r   k



  ^（5-2）

式（5-2）中 r

为第 i 个评价对象对理想对象的灰色加权关联度。

step 5、评价分析：

根据灰色加权关联度的大小，对各指标对象进行排序，可建立评价对象的 关联度，关联度越大，其评价结果越好。

二、灰色关联分析法确定最优系泊系统的设计

考虑环境阻力最大的情况下，即布放点的海水速度达到最大值 1.5m/s、风 速达到最大值 36m/s。水深在 16 米至 20 米范围内变化，编写 matlab 程序，在 1000 种设计方案中，计算得到 59 组可取的设计方案，可取的设计方案指的是 满足（5-3）条件：

0

0

0

5 16

2 h





  

 

  

（5-3）

即钢桶倾斜角度不超过 5 度，锚链末端与锚的链接处的切线方向与河床的 夹角不超过 16 度和浮标吃水深度小于 2 米。在可供选择的设计方案中，利用灰 色综合评价法综合评价浮标的吃水深度、游动区域和钢桶的倾斜角度三项指 标，运用 matlab编程，程序见附录 3，计算各种方案下，水深深度从 16 米至 20 米（水深间隔为 1 米）的关联度值并将其加总求和，排名结果见表 2：

表 2：不同海水深度下可取方案的关联度值及加总关联值 关联度值

（16m）

关联度值

（17m）

关联度值

（18m）

关联度值 (19m)

关联度值 (20m)

序号（方 案）

加总关联 度值 0.7327 0.7342 0.7328 0.7555 0.7572 4900 3.7124 0.6871 0.6901 0.689 0.6984 0.6996 4899 3.4642

0.6872 0.6823 0.6823 0.6801 0.68 4896 3.4119

0.6677 0.6701 0.6692 0.6721 0.6729 4898 3.352 0.6553 0.6609 0.6583 0.6799 0.6866 4920 3.341

0.4756 0.4824 0.4699 0.4663 0.4734 4997 2.3676

0.4659 0.4683 0.4683 0.4659 0.4757 4598 2.3441

0.4625 0.4724 0.4584 0.459 0.4677 4799 2.32

0.4583 0.466 0.4574 0.455 0.4638 4798 2.3005

0.4508 0.4518 0.4518 0.4507 0.4531 4599 2.2582 从表二最后一列可知，在序号 4900（方案）时，加总关联度值最大，则该 种设计方案效果最佳。加总关联度值等于水深分别为 16、17、18、19 和 20 米 的关联度值之和，其中序号 4900 表示某种设计方案，设计方案是：锚链的型号 为 5，长度为 20 米，重物球的质量为 3900 千克，也即最优的系泊系统的设计 方案。

5.3.2 第二个子问题的求解

不同情况下的钢桶、钢管的倾斜角度、锚链形状、浮标的吃水深度和游动 区域

当海水速度、风速固定，随着水深的增加，钢桶、钢管的倾斜角度，链环 平躺在海平面的个数和游动圆面半径在减少，浮标的吃水深度在增加，如图

（6-1）。

当海水速度、水深固定，随着风速的增大，钢桶、钢管的倾斜角度、游动 圆面半径、浮标的吃水深度在增加，链环个数基本保持不变，如图（6-2）。

当风速、水深固定，随着海水速度的增大，钢桶、钢管的倾斜角度、游动 圆面半径、浮标的吃水深度在增加,链环个数在减少,如图(6-3)。

表 6-2：不同风速下的系泊系统的参数值

风 速 （ 6

/ m s _）

风 速 （ 12

/ m s _）

风 速 （ 18

/ m s _）

风 速 （ 24

/ m s _）

风 速 （ 30

/ m s _）

风 速 （ 36

/ m s _）

钢 桶 倾 斜 角

度 0.019918 0.079643 0.179197 0.31833 0.496819 0.714047999

表 6-1：不同水深下的系泊系统的参数值

水深（16 m _）

水深（17 m _）

水深（18 m _）

水深（19 m _）

水深（20 m _）

钢桶倾斜角度 0.327074 0.321694 0.316349 0.311038 0.305762

第一节钢管倾

斜角度 0.324683 0.31935 0.314051 0.308786 0.303555

第二节钢管倾

斜角度 0.325278 0.319933 0.314622 0.309346 0.304104

第三节钢管倾

斜角度 0.325874 0.320518 0.315196 0.309908 0.304655

第四节钢管倾

斜角度 0.326473 0.321105 0.315771 0.310472 0.305207

链环平躺个数 58 53 47 42 36

游动区域半径 1.6676 1.672 1.6764 1.6808 1.6852

吃水深度 5.313989 5.210198 5.097529 4.974297 4.838082

第 一 节 钢 管

倾斜角度 0.019773 0.079063 0.17789 0.316009 0.493197 0.708844391 第 二 节 钢 管

倾斜角度 0.019809 0.079207 0.178215 0.316586 0.494098 0.71013817 第 三 节 钢 管

倾斜角度 0.019845 0.079352 0.178541 0.317165 0.495002 0.71143668 第 四 节 钢 管

倾斜角度 0.019882 0.079497 0.178869 0.317747 0.495909 0.712739948 链 环 平 躺 个

数 54 54 54 54 54 53

游 动 区 域 半

径 1.6708 1.6709 1.6709 1.6711 1.6714 1.6719

吃水深度 0.643383 1.891471 3.472208 5.183503 6.906691 8.555819087

表 6-3：不同海水速度下的系泊系统的参数值

海水速度

（0

/ m s _）

海水速度

（0.3

/ m s _）

海水速度

（0.6

/ m s _）

海水速度

（0.9

/ m s _）

海水速度

（1.2

/ m s _）

海水速度

（1.5

/ m s _）

钢 桶 倾 斜 角

度 0.322792 0.473764 0.925539 1.675626 2.719036 4.048877291 第 一 节 钢 管

倾斜角度 0.320438 0.47031 0.918799 1.663463 2.699436 4.020037259 第 二 节 钢 管

倾斜角度 0.321023 0.471169 0.920475 1.666487 2.704309 4.027208745 第 三 节 钢 管

倾斜角度 0.321611 0.472031 0.922157 1.669522 2.709201 4.03440582 第 四 节 钢 管

倾斜角度 0.3222 0.472896 0.923845 1.672569 2.71411 4.041628623 链 环 平 躺 个

数 54 54 52 47 36 19

游 动 区 域 半

径 1.6711 1.6713 1.6727 1.6766 1.6849 1.6989

吃水深度 5.232105 6.707057 9.83408 12.829214 14.939995 16.21755985

六．模型的评价和推广

6.1 模型的评价

本文就系泊系统的设计问题，建立集中质量法分析模型和灰色关联分析综 合评价模型对最优系泊系统的设计问题进行探究。

集中质量法最主要的优点是节省计算量和计算时间，有利于简化模型的求

解，但是可能会降低求解的精度，本文采用共点力的平衡条件来分析问题，将

灰色关联分析法是根据因素之间发展趋势的相似或相异程度，亦即“灰色 关联度”，作为衡量因素间关联程度的一种方法，其核心是按照一定规则确立 随时间变化的母序列，把各个评估对象随时间的变化作为子序列，求各个子序 列与母序列的相关程度，依照相关性大小得出结论。本文考虑系泊系统浮标的 吃水深度、游动区域、钢桶的倾斜角度，取各指标的最小值得到虚拟最优设 计，根据灰色加权关联度的大小确定设计的优劣之分，将定量和定性方法结合 起来，使复杂的决策问题变得简单清晰，而且计算方便，并可在一定程度上排 除决策者的主观任意性，得出比较客观的结论。

6.2 模型的改进方向

集中质量法综合分析模型的改进：

可以将钢管和链环视为刚体，考虑它们受到的重力和拉力的力矩，将系统 中物体的平衡问题视为非共点力的平衡问题。

灰色关联分析综合评价模型的改进：

可以采用线性变换的方法，重新定义关联系数。可参考石勇民《西安公路 交通大学学报》1994 年 03 期，改进的分析方法具有较高的精度，为建立模型 及确定模型参数提供了新的思路。

参考文献：

[1]司守奎，孙兆亮等，数学建模算法与应用，北京:国防工业出版社,2016.1 [2]王磊，单点系泊系统的动力学研究，中国海洋大学，P75

[3]陈宇泉，海船系泊系统的设计与安装研究，哈尔滨工程大学，U664.4

[4]童波，半潜式平台系泊系统型式及其动力特性研究，上海交通大学，U674.381 [5]郭小天 ，漂 浮式 潮流电站 弹性 系泊系 统设计研 究， 哈尔滨 工程大学 ， TV744;U653.2

[6]何静， _FPSO 悬式锚腿系泊系统的锚系设计研究，武汉理工大学，U653.2 [7]石勇民,西安公路交通大学学报

http://xueshu.baidu.com/s?wd=paperuri:(4e18149703a5ba5ab92a5c68060ce8

d6)&filter=sc_long_sign&tn=SE_baiduxueshuwise_y5e9t2rx&sc_us=11893325

829204560143

附录 附录一

（matlab 编写的求解问题一、二、三的代码）

%运行程序的代码，调用各个自己编写的力函数，可调整 Problem1、2、3、4 的 true 或 false 来针对不同问题进行求解 clc,clear

SeaDensity = 1025; %海底密度 g = 9.8; %重力加速度 SeaHeghit = 18; %海底深度

Buoy_r = 1; %浮标的半径 Buoy_h = 2; %浮标的高 Buoy_m = 1000; %浮标的质量

GangGuan_r = 0.025; %钢管的半径 GangGuan_h = 1; %钢管的高 GangGuan_m = 10; %钢管的质量

GangGuan_Theta = zeros(4,1); %初始化4个钢管的角度 GangGuan_LaLi = zeros(4,1); %初始化4个钢管的拉力

GangTong_r = 0.15; %钢桶的半径 GangTong_h = 1; %钢桶的高 GangTong_m = 100; %钢桶和设备的质量

ZhongWuQiu_m = 1200; %重物球的质量(问题一）

%ZhongWuQiu_MiDu = 7900; %钢铁的密度

LianHuanH = [0.78, 0.105, 0.12, 0.15, 0.18]; %不同型号锚链的长度 LianHuan_P = [3.2, 7, 12.5, 19.5, 28.12]; %锚链的单位长度的质量

LianHuan_Theta = zeros(300,1); %初始化链环的角度 LianHuan_LaLi = zeros(300,1); %初始化链环的拉力

Theta = zeros(300, 1); %初始化用来保存各个器件的角度 Num = 1; %用来记录目前保存到哪个角度值

Problem1 = false; %问题1求解 Problem2 = false; %问题2求解

Problem3 = false; %问题3第一小问求解（确定系泊系统的设计）

Problem4 = true; %问题4第二小问求解（分析不同情况下各个值）

global G_ShuiLiuLi; %设置水流力为全局变量，提供给F_ShuiLiuLi函数用

G_FengSu = 36;

G_ShuiLiuLi = 0;

%初始化风速的条件 FengSu_ChuShiZhi = 12;

FengSu_ZhongZhiZhi = 12;

%初始化水流力的条件 ShuiLiuLi_ChuShiZhi = 0;

ShuiLiuLi_ZhongZhiZhi = 0;

%初始化海底深度的条件

SeaHeghit_ChuShiZhi = 18; %海底深度 SeaHeghit_ZhongZhiZhi = 18; %海底深度

%初始化链环的条件

LianHuan_XH_ChuShiZhi = 2;

LianHuan_XH_ZhongZhiZhi = 2;

%初始化重物球的条件

ZhongWuQiu_m_ChuShiZhi = 1200;

ZhongWuQiu_m_JingDu = 1;

ZhongWuQiu_m_ZhongZhiZhi = 1200;

ZhongWuQiu_TiaoJie = true;

%初始化链环的条件

LianHuan_Len_ChuShiZhi = 22.05;

LianHuan_Len_JingDu = 1;

LianHuan_Len_ZhongZhiZhi = 22.05;

LianHuan_TiaoJie = false;

if (Problem1 == true)

Buoy_hx_JingDu = 0.00001; %对于问题1和问题3要求浮标海底精度高使得计算的海底深度接近于18

Buoy_hx_ChuShiZhi = 0.73; %在各个问题1的限制条件下，浮标海底的初始值不会低于0.7优化运算过程，问题3时可 适当调整降低运算次数

elseif (Problem2 == true)

Buoy_hx_JingDu = 0.001; %对于问题2浮标海底精度不要求那么高，减少可接受精度提高运算速度 Buoy_hx_ChuShiZhi = 0.9;

ZhongWuQiu_m_ChuShiZhi = 1200;

ZhongWuQiu_m_JingDu = 100;

ZhongWuQiu_m_ZhongZhiZhi = 10000;

elseif (Problem3 == true)%调节链环的总长度

Buoy_hx_JingDu = 0.01; %降低精度提高运算速度 Buoy_hx_ChuShiZhi = 0.6;

LianHuan_Len_ChuShiZhi = 16;

LianHuan_Len_JingDu = 1;

LianHuan_Len_ZhongZhiZhi = 25;

ZhongWuQiu_m_ChuShiZhi = 2000;

ZhongWuQiu_m_JingDu = 100;

ZhongWuQiu_m_ZhongZhiZhi = 3900;

LianHuan_XH_ChuShiZhi = 1;

LianHuan_XH_ZhongZhiZhi = 5;

SeaHeghit_ChuShiZhi = 16;

SeaHeghit_ZhongZhiZhi = 20;

Answer = zeros(5000, 5); %初始化最终5000个结果 elseif (Problem4 == true)

Buoy_hx_JingDu = 0.0001;

Buoy_hx_ChuShiZhi = 1.6;

SeaHeghit = 16.79310718;

G_FengSu = 24.16760985;

G_ShuiLiuLi = 0.853493409;

Answer = zeros(5, 8);

FengSu_ChuShiZhi = 24.16760985;

FengSu_ZhongZhiZhi = 24.16760985;

SeaHeghit_ChuShiZhi = 16.79310718;

SeaHeghit_ZhongZhiZhi = 16.79310718;

ShuiLiuLi_ChuShiZhi = 0;

ShuiLiuLi_ZhongZhiZhi = 1.5;

%设置最优系泊系统参数 LianHuan_XH_ChuShiZhi = 5;

LianHuan_XH_ZhongZhiZhi = 5;

LianHuan_Len_ChuShiZhi = 20;

LianHuan_Len_JingDu = 1;

LianHuan_Len_ZhongZhiZhi = 20;

ZhongWuQiu_m_ChuShiZhi = 3900;

ZhongWuQiu_m_JingDu = 1;

end;

for G_ShuiLiuLi = ShuiLiuLi_ChuShiZhi:0.3:ShuiLiuLi_ZhongZhiZhi;

for G_FengSu = FengSu_ChuShiZhi:6:FengSu_ZhongZhiZhi;

for SeaHeghit = SeaHeghit_ChuShiZhi:SeaHeghit_ZhongZhiZhi; %海底深度

for LianHuan_XH = LianHuan_XH_ChuShiZhi:LianHuan_XH_ZhongZhiZhi; %锚链的型号 LianHuan_h = LianHuanH(LianHuan_XH); %链环的单位长度

%LianHuan_r = 0.0095; %链环的半径

LianHuan_m = LianHuan_h * LianHuan_P(LianHuan_XH); %链环的质量

for LianHuan_Len = LianHuan_Len_ChuShiZhi:LianHuan_Len_JingDu:LianHuan_Len_ZhongZhiZhi; %链环 的总长度

LianHuan_Num = round(LianHuan_Len / LianHuan_h); %计算链环的数量(四舍五入取整）

for ZhongWuQiu_m = ZhongWuQiu_m_ChuShiZhi:ZhongWuQiu_m_JingDu:ZhongWuQiu_m_ZhongZhiZhi; %第二题_调 整重物球的质量让角度符合题目范围

%通过调整Buoy_hx让算的得高度接近于海底深度

for Buoy_hx = Buoy_hx_ChuShiZhi:Buoy_hx_JingDu:2; %浮标在海底的高度 if (Problem3 == true || Problem4 == true)

Buoy_ShuiLiuLi = F_ShuiLiuLi(Buoy_r, Buoy_h, Buoy_hx, 0); %浮标的水流力 else

Buoy_ShuiLiuLi = 0;

end; %if (Problem3 == true)

Hsum = Buoy_hx; %实际计算的海底高度 Hxsum = 0; %游动区域

Buoy_FuLi = F_FuLi(Buoy_r, Buoy_hx); %浮标的浮力 Buoy_FengLi = F_FengLi(Buoy_r, Buoy_h, Buoy_hx); %浮标的风力

Buoy_LaLi = sqrt((Buoy_FuLi - Buoy_m * g)^2 + (Buoy_FengLi + Buoy_ShuiLiuLi)^2); %浮标的拉 力

Buoy_Theta = atan((Buoy_FengLi + Buoy_ShuiLiuLi)/(Buoy_FuLi - Buoy_m * g)) + pi; %浮标的 theta角的值

if (Buoy_FuLi < Buoy_m * g) %如果浮力<重力，则浮标下沉，不考虑 break;

end; %if (Buoy_Theta < pi)

Num = 1; %初始化Num的值

Theta(Num,1) = Change_Theta(Buoy_Theta); %上一个夹角的值 LaLi = Buoy_LaLi; %上一个拉力的值

GangGuan_FuLi = F_FuLi(GangGuan_r, GangGuan_h); %钢管的浮力 if (Problem3 == true || Problem4 == true)

GangGuan_ShuiLiuLi = F_ShuiLiuLi(GangGuan_r, GangGuan_h, 0, Theta(Num,1)); %钢管的水流

力

else

GangGuan_ShuiLiuLi = 0;

end; %if (Problem3 == true)

for i = 1:4;

[GangGuan_Theta(i,1)] = F_theta(LaLi, Theta(Num,1), GangGuan_FuLi, GangGuan_m * g, GangGuan_ShuiLiuLi); %钢管的角度

[GangGuan_LaLi(i,1)] = F_LaLi(GangGuan_FuLi, LaLi, Theta(Num,1), GangGuan_m * g, GangGuan_ShuiLiuLi); %钢管的拉力

Num = Num + 1;

Theta(Num,1) = Change_Theta(GangGuan_Theta(i,1)); %上一个夹角的值 LaLi = GangGuan_LaLi(i,1); %上一个拉力的值

Hsum = Hsum + cos(Theta(Num - 1,1)) * GangGuan_h; %计算四个钢管竖直方向的高度 Hxsum = Hxsum + sin(Theta(Num - 1,1)) * GangGuan_h; %计算四个钢管水平方向的长度 end; %for i = 1:4;

if (Problem3 == true || Problem4 == true)

GangTong_ShuiLiuLi = F_ShuiLiuLi(GangTong_r, GangTong_h, 0, Theta(Num,1)); %钢桶的水流 力

else

GangTong_ShuiLiuLi = 0;

end; %if (Problem3 == true)

GangTong_FuLi = F_FuLi(GangTong_r, GangTong_h); %钢桶的浮力

%GangTong_G = GangTong_m * g + ZhongWuQiu_m * g - SeaDensity * g * ZhongWuQiu_m / ZhongWuQiu_MiDu; %钢桶的重力和重物球的重力和浮力

GangTong_G = GangTong_m * g + ZhongWuQiu_m * g; %钢桶的重力和重物球的重力 GangTong_Theta = F_theta(LaLi, Theta(Num,1), GangTong_FuLi, GangTong_G,

GangTong_ShuiLiuLi); %钢桶的角度

GangTong_LaLi = F_LaLi(GangTong_FuLi, LaLi, Theta(Num,1), GangTong_G, GangTong_ShuiLiuLi); %钢桶的拉力

Num = Num + 1;

Theta(Num,1) = Change_Theta(GangTong_Theta); %上一个夹角的值 LaLi = GangTong_LaLi; %上一个拉力的值

Hsum = Hsum + cos(Theta(Num - 1,1)) * GangTong_h; %计算钢桶竖直方向的高度 Hxsum = Hxsum + sin(Theta(Num - 1,1)) * GangTong_h; %计算钢桶水平方向的长度

%LianHuan_FuLi = F_FuLi(LianHuan_r, LianHuan_h, 0); %链环的浮力 LianHuan_FuLi = 0; %不考虑链环的浮力

LianHuan_ShuiLiuLi = 0; %不考虑链环的水流力

[LianHuan_Theta(i,1)] = F_theta(LaLi, Theta(Num,1), LianHuan_FuLi, LianHuan_m * g, LianHuan_ShuiLiuLi);%链环的角度

[LianHuan_LaLi(i,1)] = F_LaLi(LianHuan_FuLi, LaLi, Theta(Num,1), LianHuan_m * g, LianHuan_ShuiLiuLi); %链环的拉力

Num = Num + 1;

Theta(Num,1) = Change_Theta(LianHuan_Theta(i,1)); %上一个夹角的值 LaLi = LianHuan_LaLi(i,1); %上一个拉力的值

Hsum = Hsum + cos(Theta(Num - 1,1)) * LianHuan_h; %计算链环竖直方向的高度 Hxsum = Hxsum + sin(Theta(Num - 1,1)) * LianHuan_h; %计算链环水平方向的长度 end; %for i = 1:LianHuan_Num;

if (Hsum > SeaHeghit) %如果计算出的高度小于等于海底实际高度则跳出 break;

end;%if (Hsum <= SeaHeghit) end; %for_浮标在海面外的高度

if (Problem1 == true) %问题一下显示所要求的的数值 for i = 1:300 %显示各个器件拉力的角度 if (Theta(i, 1) == 0)

break;

end; %if (Theta(i) == 0) end; %for i = 1:300

Theta(1:i-1, 1)

Hsum %显示实际计算的海底高度 ChiShuiShenDu = Buoy_hx %显示浮标在海底的高度

MaoLianJiaJiao = 90 - (Theta(LianHuan_Num + 6) * 180 / pi) %显示锚链末端与锚的链接处的 切线方向与海床的夹角

GangTongQingXieJiao = Theta(5) * 180 / pi %显示钢桶的倾斜角度（钢桶与竖直线的夹角）

elseif (Problem2 == true) %问题二下显示所要求的的数值 if (Buoy_hx >= 2) %如果浮标完全沉入水底则停止运算 break;

end; %if (Buoy_hx > 2)

if (90 - (Theta(LianHuan_Num + 6) * 180 / pi) > 16 || Theta(5) * 180 / pi > 5) continue;

end; %if (90 - (Theta(LianHuan_Num + 6) * 180 / pi) > 16 || Theta(5) * 180 / pi > 5)

ZhongWuQiu_m %显示重物球的质量 ChiShuiShenDu = Buoy_hx %显示浮标在海底的高度

MaoLianJiaJiao = 90 - (Theta(LianHuan_Num + 6) * 180 / pi) %显示锚链末端与锚的链接处的 切线方向与海床的夹角

GangTongQingXieJiao = Theta(5) * 180 / pi %显示钢桶的倾斜角度（钢桶与竖直线的夹角）

elseif (Problem3 == true) i = SeaHeghit - 16;

j = LianHuan_XH - 1;

k = LianHuan_Len - 16;

l = (ZhongWuQiu_m - 1900) / 100;

AnswerNum = i * 5 * 10 * 20 + j * 10 * 20 + k * 20 + l;

Answer(AnswerNum, 1) = Buoy_hx; %保存吃水深度

MaoLianJiaJiao = 90 - (Theta(LianHuan_Num + 6) * 180 / pi); %保存锚链末端与锚的链接处 的切线方向与海床的夹角

GangTongQingXieJiao = Theta(5) * 180 / pi; %保存钢桶的倾斜角度（钢桶与竖直线的夹 角）

Answer(AnswerNum, 3) = GangTongQingXieJiao; %保存钢桶倾斜角度 Answer(AnswerNum, 2) = Hxsum; %保存游动区域

if (GangTongQingXieJiao > 5 || MaoLianJiaJiao > 16) %确定是否满足系统 Answer(AnswerNum, 4) = 0;

else

Answer(AnswerNum, 4) = 1;

end;

Answer(AnswerNum, 5) = AnswerNum; %确定答案的序号 elseif (Problem4 == true)

%i = SeaHeghit - 15;

%i = G_FengSu / 6;

i = G_ShuiLiuLi * 10 / 3 + 1;

Answer(i, 1) = Theta(5) * 180 / pi; %记录钢桶 Answer(i, 2) = Theta(1) * 180 / pi; %记录钢管 Answer(i, 3) = Theta(2) * 180 / pi;

Answer(i, 4) = Theta(3) * 180 / pi;

Answer(i, 5) = Theta(4) * 180 / pi;

%计算沉在海底的链环个数 HaiDiNum = 0;

HaiDiShenDu = 0;

XNum = Num - 1;

for j = 1:LianHuan_Num;

HaiDiShenDu = HaiDiShenDu + LianHuan_h * cos(Theta(XNum,1));

XNum = XNum - 1;

if (HaiDiShenDu <= 0) HaiDiNum = j;

elseif (abs(HaiDiNum) < 1e-8) break;

end;

end;

Answer(i, 6) = HaiDiNum;

Answer(i, 7) = Buoy_hx;

Answer(i, 8) = Hxsum;

end; %if (Problem1 == true)

if (Problem3 == false) %问题三不做图 if (Num <= 1)

continue;

end;

Num = Num - 1; %因为最后一个角度是算的锚的角度 LianHuan_x = zeros(LianHuan_Num + 1, 1);

LianHuan_y = zeros(LianHuan_Num + 1, 1);

LianHuan_yy = zeros(LianHuan_Num + 1, 1); %画出实际y轴高度 HaiDiNum = 0;

for i = 2:LianHuan_Num + 1

LianHuan_x(i) = LianHuan_x(i-1) + LianHuan_h * sin(Theta(Num,1));

LianHuan_y(i) = LianHuan_y(i-1) + LianHuan_h * cos(Theta(Num,1));

if (LianHuan_y(i) > 0)

LianHuan_yy(i) = LianHuan_y(i);

else

LianHuan_yy(i) = 0;

HaiDiNum = i;

end;

Num = Num - 1;

end;

hold on;

plot(LianHuan_x', LianHuan_yy', '-');

axis('equal');

%绘制钢桶图形

GangTong_x = zeros(2,1);

GangTong_y = zeros(2,1);

GangTong_x(1,1) = LianHuan_x(LianHuan_Num + 1,1);

GangTong_y(1,1) = LianHuan_y(LianHuan_Num + 1,1);

GangTong_x(2,1) = GangTong_x(1,1) + GangTong_h * sin(Theta(Num,1));

GangTong_y(2,1) = GangTong_y(1,1) + GangTong_h * cos(Theta(Num,1));

Num = Num - 1;

plot(GangTong_x, GangTong_y, '*-', 'LineWidth', 3)

%绘制钢管图形

GangGuan_x = zeros(5,1);

GangGuan_y = zeros(5,1);

GangGuan_x(1,1) = GangTong_x(2,1);

GangGuan_y(1,1) = GangTong_y(2,1);

for i = 2:5

GangGuan_x(i,1) = GangGuan_x(i-1,1) + GangGuan_h * sin(Theta(Num,1));

GangGuan_y(i,1) = GangGuan_y(i-1,1) + GangGuan_h * cos(Theta(Num,1));

Num = Num - 1;

end;

plot(GangGuan_x, GangGuan_y, '*-', 'LineWidth', 2)

%绘制浮标图形

FuBiao_x = zeros(2,1);

FuBiao_y = zeros(2,1);

FuBiao_x(1,1) = GangGuan_x(5,1);

FuBiao_y(1,1) = GangGuan_y(5,1);

FuBiao_x(2,1) = FuBiao_x(1,1);

FuBiao_y(2,1) = FuBiao_y(1,1) + Buoy_hx;

plot(FuBiao_x, FuBiao_y, '-', 'LineWidth', 20)

title('不同风速下系泊系统的形状');

legend('锚链','钢桶','钢管','浮标', 0)

hold off;

if (Problem1 == true || Problem3 == true)

YouDongQuYu = FuBiao_x(2,1) %浮标的游动区域，浮标的平衡点及浮标的x轴坐标 end; %if (Problem1 == true || Problem3 == true)

end;%if (Problem3 == false)

end; %for ZhongWuQiu_m = ZhongWuQiu_m_ChuShiZhi:ZhongWuQiu_m_JingDu:ZhongWuQiu_m_ZhongZhiZhi;

end;%for LianHuan_Len = LianHuan_Len_ChuShiZhi:LianHuan_Len_JingDu:LianHuan_Len_ZhongZhiZhi;

end;%for LianHuan_XH = LianHuan_XH_ChuShiZhi:LianHuan_XH_ZhongZhiZhi;

end; %for SeaHeghit_ChuShiZhi:SeaHeghit_ZhongZhiZhi;

end;%for G_FengSu = FengSu_ChuShiZhi:6:FengSu_ZhongZhiZhi;

end;%for G_ShuiLiuLi = ShuiLiuLi_ChuShiZhi:6:ShuiLiuLi_ZhongZhiZhi;

if (Problem4 == true) Answer'

end; %if (Problem4 == true) if (Problem3 == true)

xlswrite('text', Answer); %将最终结果导出到excel文件中，

end; %if (Problem3 == true) %筛选数据时，需再将数据完整的复制到in.txt中，运行数据筛选.exe，再将 out.txt复制到任意excel文件中即可

YunXingJieGuo = '运行完成'

%风力函数,保存在 F_FengLi.m 中

function [y] = F_FengLi(r, h, hx)

%计算风力，r 为物体的半径，h 为物体的高，hx 为物体在海底的高度

global G_FengSu

FengSu = G_FengSu; %风速的大小

y = 0.625 * (h - hx) * 2 * r * FengSu^2; %F = 0.625 * S * v"2(N)

end

function [y] = F_FuLi(r, hx)

%计算圆柱体浮力，r 为物体的半径，hx 为物体在海底的高

SeaDensity = 1025; %海底密度

g = 9.8; %重力加速度

y = SeaDensity * g * (pi * r^2 * hx); %F = pgv;

end

%拉力函数,保存在 F_LaLi.m 中

function [y] = F_LaLi(F_FuLi, F_LaLi, theta, G, F_ShuiLiuLi)

%计算该物体对下一个物体的拉力，F_FuLi 为该物体的浮力，F_LaLi 为上一个物体对该物体的拉力

%theta 为上一个物体对该物体的夹角，G 为该物体的重力，F_ShuiLiuLi 为该物体的水流力

y = sqrt((F_FuLi + F_LaLi * cos(theta) - G)^2 + (F_LaLi * sin(theta) + F_ShuiLiuLi)^2);%F_LaLi = sqrt(Ty"2 + Tx^2) Tx 和 Ty 分别

是水平和竖直方向的合力

end %function [y] = F_LaLi(F_FuLi, F_LaLi, theta, G, F_ShuiLiuLi)

%水流力函数,保存在 F_ShuiLiuLi.m 中

function [y] = F_ShuiLiuLi(r, h, hx, theta)

%计算水流力，r 为物体的半径，h 为物体的高，hx 为物体在海底的高度

global G_ShuiLiuLi

ShuiLiuSu = G_ShuiLiuLi; %水流力的大小

y = 374 * hx * 2 * r * ShuiLiuSu^2 * cos(theta); %F = 374 * S * v"2

end

%计算角度的函数,保存在 F _Theta.m 中

function [y] = F_theta(F_LaLi, theta, F_FuLi, G, F_ShuiLiuLi)

%计算该物体对下一个物体的拉力方向角，F_LaLi 为上一个物体对该物体的拉力，theta 为上一个物体对该物体的夹角，

%F_FuLi 为该物体的浮力，G 为该物体的重力，F_ShuiLiuLi 为该物体的水流力

Tx = F_LaLi * sin(theta) + F_ShuiLiuLi; %水平方向的合力

Ty = F_FuLi + F_LaLi * cos(theta) - G; %竖直方向的合力

if (Tx > 0 && abs(Ty) < 1e-8) %考虑竖直方向的合力为 0 的情况

y = (3/2) * pi;

elseif (Tx < 0 && abs(Ty) < 1e-8)

y = (1/2) * pi;

end; %if (Tx > 0 && abs(Ty) < 1e-8)

if (abs(Ty) >= 1e-8) %考虑水平方向的合力或竖直方向的合力为负数的情况

y = atan( Tx/Ty );

if (Tx >= 0 && Ty > 0)

y = y + pi;

elseif (Tx > 0 && Ty < 0)

y = y + 2 * pi;

elseif (Tx < 0 && Ty > 0)

y = y + pi;

elseif (abs(Tx) < 1e-8 && Ty < 0)

y = 0;

end; %if (Tx >= 0 && Ty > 0)

end; %if (abs(Ty) >= 1e-8)

end %function [y] = F_theta(F_LaLi, theta, F_FuLi, G, F_ShuiLiuLi)

%更改角度函数,保存在 Change_Theta.m 中

function [y] = Change_Theta(Theta)

%上一个物体对该物体拉力的倾斜角度变换成该物体受到上一物体对它的倾斜角度（转换受力分析的对象）

if (Theta <= pi)

y = Theta + pi;

else

y = Theta - pi;

end; %if (Theta <= pi)

end %function [y] = Change_Theta(Theta)

附录二

（matlab 实现灰色关联分析模型）

a=[0.74423 0.75928 0.78325 0.81611 0.85879

18.48150818 18.31874982 17.41479432 16.00962318 14.46753747

6.234253012 6.023759016 5.716482983 5.343101259 4.925652369];%在风速 24m/s,水流速度 1m/s,水深 18m，锚链 22.05m,重物球质量为 1200kg 条件下，矩阵 a 第 j 列表示型号为 j 的锚链，浮标的吃水深度、游动区域、钢桶倾斜角度。

for i=1:3%成本型指标标准化

a(i,:)=(max(a(i,:))-a(i,:))/(max(a(i,:))-min(a(i,:)));

end

[m,n]=size(a);

cankao=max(a')'%求参考数列的取值

t=repmat(cankao,[1,n])-a;%求参考数列与每一个序列的差 mmin=min(min(t));%计算最小差

mmax=max(max(t));%计算最大差 rh0=0.5;%分辨系数

xishu=(mmin+rh0*mmax)./(t+rh0*mmax)%计算灰色关联系数 guanliandu=mean(xishu)%取等系数权重，计算关联度

[gsort,ind]=sort(guanlaindu,'descend')%对关联度降序排列