当前位置：文档库 › 2016数学建模A题系泊系统设计

2016数学建模A题系泊系统设计

系泊系统的设计

摘要

对于问题一，建立模型一，已知题目给出的锚链长度与其单位长度的质量，得到悬链共210环。对各节锚链，钢桶，四节钢管受力分析得出静力平衡方程，使用分段外推法，可以得到静力平衡下的迭代方程。其中锚对锚链的拉力大小方向为输入变量，迭代的输出变量为浮标的位置和对钢管的拉力，在给定的风速下，输入和输出满足关系2)2(25.1cos 水v h T -=α，αθcos cos 11T T =，通过多层搜索算法得出最符合的输入输出值，即可得到给定风速下浮标的吃水深度，浮标拉力、锚链与海床夹角。利用MATLAB 软件编程求解模型得到：风力12m/s 时，钢桶与竖直方向上的角度1.9863度，从下往上四节钢管与竖直方向夹角为1.9652度、1.9592度、1.9532度、1.9472度，浮标吃水0.7173m ，以锚为圆心浮标的游动区域16.5125m ，锚链末端切线与海床的夹角3.8268度。风力24m/s 时，锚链形状，钢桶与竖直方向上的夹角3.9835度，从下往上四节钢管与竖直方向夹角为3.9420度、3.9301度、3.9183度、3.9066度，浮标吃水0.7244m ，以锚为圆心浮标的游动区域18.3175m 。锚链末端切线与海床夹角15.9175度。

对于问题二的第一小问，使用模型一求解，当风速36m/s 时，锚链末端切线与海床夹角26.3339度，浮标吃水0.7482m ，浮标游动区域为以锚为圆心半径为18.9578m 的圆形区域，从下往上四节钢管与竖直方向倾斜角度为8.4463度、8.4225度、8.3989度、8.3753度，钢桶与竖直方向倾斜角度为8.5294度。为满足问题二的要求，在模型一的基础上把重物球质量作为变量，建立模型二，将钢桶倾斜角小于5度和锚链前端夹角小于16度当做两个约束条件，通过MATLAB 编程求解得到满足约束条件要求的重物球质量取值范围为3700kg 到5320kg 。

对于问题三，首先取不同水深、水速、风速三种情况，建立模型三，即在模型一的基础上增加水流对系统产生的影响。在三种情况下，找到合适的锚链型号、锚链长度，重物球质量，对吃水深度、游动区域、钢桶的倾斜角三个目标进行优化达到最小。通过MATLAB 编程实现该模型三得到结果：选用Ⅲ型锚链，锚链长度为27.24m ，重物球质量为2580kg 。

关键词：平面静力系分析多层搜索算法遗传算法逐步外推法多目标优化

一、问题重述

近浅海观测网的传输节点由三部分组成：浮标系统、系泊系统和水声通讯系统，如下图所示：浮标系统可简化为一个圆柱体，质量为1000kg 。系泊系统由五部分组成：钢管、钢桶、重物球、电焊锚链和特制的抗拖移锚。锚的质量为600kg ，常用型号及其参数在附表中列出。钢管共4节，每节质量为10kg 。锚链末端与锚的链接处的切线方向与海床的夹角不超过?16。水声通讯系统安装在圆柱形钢桶内，设备和钢桶总质量为100kg 。钢桶上接第4节钢管，下接电焊锚链，钢桶的倾斜角度不超过?5。为控制钢桶的倾斜角度，钢桶与锚链链接处可悬挂重物球。

传输节点示意图

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

1.某型传输节点选用锚链型号II ，长度2

2.05m ，选用的重物球的质量为1200kg 。现将该型传输节点布放在水深18m 、海床平坦、海水密度为

33m /kg 10025.1?的海域。假设海水静止，分别计算海面风速为12m/s 和24m/s 时钢桶和各节钢管的倾斜角度、锚链形状、浮标的吃水深度和游动区域。

2.在问题1的假设下，计算海面风速为36m/s 时钢桶和各节钢管的倾斜角度、锚链形状和浮标的游动区域。请调节重物球的质量，使得钢桶的倾斜角度不超过?5，锚链在锚点与海床的夹角不超过?16。

3.由于潮汐等因素的影响，布放海域的实测水深介于20m -16m 之间。布放点的海水速度最大为1.5m/s 、风速最大为36m/s 。请给出考虑风力、水流力和水深情况下的系泊系统设计，分析不同情况下钢桶、钢管的倾斜角度、锚链形状、浮标的吃水深度和游动区域。

注：近海风荷载的近似公式为)(v 625.02N S F ?=，S 为物体在风向法平面的投影面积)m (2，v 为风速)(m/s 。近海水流力的近似公式为)(v 3742N S F ?=，其中S 为物体在水流速度法平面的投影面积)m (2，v 为水流速度)(m/s 。

表注：长度是指每节链环的长度。

二、问题分析

对于问题一，所选取锚链的长度和型号是给定的，所以其锚链的节数、每一节的长度、每一节的质量都是已知的，该问中选用II型的锚链，根据题目中的附表求得锚链分为210节，每节重0.735kg。根据题目所给图示可知该系泊系统的构成是：在海底的锚连接锚链的下端，在锚链的上端也就是钢桶的下端节点处挂有重物球，钢桶的上端连接4个钢管，钢管的上端连接浮标。可以假设链环之间、锚链和其他结构的连接处、钢管之间的连接为铰接，且所有结构都是质量分布均匀的，故每个链环和钢桶和钢管可以抽象为体积不同质量不同的铰接圆柱体，在独立分析出每个圆柱体所受浮力来之后，可以再把各圆柱体抽象为有质量的铰接直杆。对与锚连接的第一个链环进行受力分析如图1、图2、图3、图4所示：

图1 图2

图3 图4

对于该链环，假设和锚连接的节点的力与该链环方向一致，需要满足力平衡与力矩平衡。列出所有杆的平衡方程，可以通过方程迭代求解出在不同拉力1T作

用时的对应的浮标的位置，利用浮标位置和风力、拉力的关系式，可以再通过优化算法得到满足关系式的最优解分别求解风速为12m/s和24m/s的情况下的。

对于问题二，可以继续利用在问题一中建立的模型求解出风速在36m/s的情况下各个需要求的量。接着在问题一模型的基础上进行修改，将重物球的质量改为变量当做模型二，把该问题中对钢桶倾斜角小于5度和锚链前端夹角小于16度的要求当做两个约束条件来限制重物球质量的范围，把在充分大范围内按一定精度取值的重物球质量带入模型二中，得到一系列对应的钢桶倾斜角和锚链前端夹角值，取两者都符合约束条件的值对应的一组质量值，该组值的最小值和最大值即为可调节重物球质量的下限和上限。

对于问题三，由于海水的速度可能是由于风力的影响产生的，所以可以认为海水速度和风速是相同方向的。考虑到水里的部分都需要受到水流力的影响所以在受力分析时需要添加水平方向的作用力，如图5、图6所示：

图5 图6

由于锚链的体积忽略，故锚链的部分没有水流力的影响，对于其他部分可写出每个杆的力平衡方程，力矩平衡方程。由于水深、水流速、风速都是在一定范围内的变量，所以要寻找合适的锚链型号和链长度以及合适的重物球质量，来使得吃水深度、钢桶倾斜角、和游动区域尽可能小。该问题属于多目标优化问题，其中优化的目标为吃水深度、钢桶倾斜角、游动区域求最小；变量有水深、水流速、风速、锚链型号、链长度，其中的水深、水流速、锚链要分别取不同的情况，每个优化目标在不同情况下得到值的和最小所对应的设计方案即为最优方案。

三、模型假设

1.假设重物球的体积可以忽略

2.假设浮标系统在水中的倾斜可以忽略

3.假设锚链是不可压缩或伸长的

4.假设锚链在水中的体积可以忽略

5.假设各连接点都是铰连接方式

6.假设忽略波浪力的作用

四、符号说明

v水的流速

水

v风速

风

ρ海水密度

海

h浮标自身吃水深度

h海深度

海

l钢管长度

管

l钢桶长度

桶

V钢管体积

管

V钢桶体积

桶

h浮标高度

浮

V浮标体积

浮

h浮标吃水

m钢桶质量

桶

m钢管质量

管

m重物球质量

重

m每节锚链质量

x第n个节点横坐标

y第n个节点纵坐标

g重力加速度

l单节链的长度

链

T锚对链的拉力

θ锚链末端与海床的夹角

T浮标对链的拉力

α浮标拉力和海平面夹角

θ钢桶与竖直线的夹角

桶

x最大游动区域半径

max

五、模型建立与求解

5.1问题一

5.1.1模型一的建立

首先对锚链进行受力分析，如图7所示

图7

拉力1T 是沿杆方向的拉力，重力作用于杆的中心产生力矩，拉力2T 需要平衡拉力1T 和重力产生的力矩，才能使所有杆满足力平衡和力矩平衡：

0,0==∑∑yn xn F F 0,0==∑∑yn xn

M M

其中n 为所有杆个数，这里根据题目里提供的数据，其锚链的节数、每一节的长度、每一节的质量都是已知的，该问中选用II 型的锚链，根据题目中的附表求得锚链分为210节，整个抽象杆除此之外还有一节钢桶和四节钢杆组成。因此215=n 。

对于前209个杆满足(5.2)式的迭代关系：

1111111111sin cos 00θθθθ=====T T T T y x y x （5.1）

)

1()

1(1)1()1(arctan 21

+++++=+==n x n y n l yn n y xn

n x T T

g m T T T T θ （5.2）

其中（5.1）式作为初始的条件。

在接入重物球段也就是210个杆的迭代表达式不同于（5.2）为：

q l yn n y m g m T T ++=+2

)1(

钢桶的为：

g m T T t yn n y 2

)1(+=+

钢管的为：

g m T T g yn n y 2

)1(+=+

以上各部分都满足（5.3）式

n n n n n l y y l x x θθsin cos 11+=+=++ （5.3）

通过（5.3）迭代求解出每个1θ和1T 对应的最终的节点位置即浮标下底面中

心的坐标),(216216y x 和该点的作用力),(216216y x T T ，可以将以上迭代过程省略，写

成

),(],,,[11216216216216T f T T y x y x θ=

（5.4）

对整个铰接杆结构受力分析可以得到水平方向上不受外力作用，故由于受力平衡得到：

αθcos cos 11T T = （5.5）

在忽略浮标倾斜的情况下，浮标对杆的拉力的水平分量大小等于浮标受到的风力

)2(25.1cos 水v h T -=α （5.6）

其中

T h h s 海πρα

sin +

= （5.7）

h y -=18216 （5.8）

这里需要完成：求一组),(11T θ可以满足（5.5）、（5.6）、（5.7）、（5.8）四式，联立可得到三个方程

θαπραθcos cos )2025.1(cos )sin 2(25.1cos 112

2162

11T T v y T v g

T h T s =-=--=风

风

海（5.9）

这里需要完成：求一组),(11T θ可以满足（5.9）式。利用优化的思想该问题

可以转化为

???

???????????-+--+---αθαπραθcos cos )2025.1(cos )sin 2(25.1cos min 11221621T T v y T v g T h T s 风风海（5.10）

有约束

??≤≤≤≤10000090

0..11T t s θ （5.11）

（5.10）中各个量可以通过下面式子求得

),(],,,[11216216216216T f T T y x y x θ=

5.1.2模型一的求解

对于模型一所建立的优化模型，有多种求解方法，考虑到各种优化算法的效率以及准确度。这里采用多层搜索算法求解。

多层搜索算法的主要步骤有如下几步：

1.确定精度系数和初始搜索空间

2.将搜索空间按精度系数分块分出多个点，找出取值最优的点

3.在最优点的左右取精度系数分块大小作为新的搜索空间

4.重复2到4直到满足条件

按照此算法步骤，使用MATLAB 进行编程实现，代码见附录。运行程序可以得到结果如下：

风速12m/s 时 3.82681=θ，浮标吃水0.7173米，游动区域为以海面上锚纵坐标为圆心半径为16.5125米的圆形区域，从下往上四节钢管在竖直方向倾斜角度为1.9652度、1.9592度、1.9532度、1.9472度。钢桶在竖直方向倾斜角度为1.9863度。锚链形状如图8所示

图8

风速24m/s 时15.91751=θ，吃水0.7244m ，游动区域为以海面上锚纵坐标为圆心半径为18.3175的圆形区域，从下往上四节钢管在竖直方向倾斜角度为3.9420度、3.9301度、3.9183度、3.9066度钢桶在竖直方向倾斜角度为3.9835度。锚链形状如图9所示

图9

问题一中12m/s的系泊系统简单图示如图10：

图10

5.2问题二

5.2.1模型二的建立

对于第一小问可以直接使用模型一进行求解。对于求重物球的质量范围，需要重物球为变量所以可以根据所建立的模型一，改进得到模型二。在模型一中的重物球是固定的，在模型二中重物球的质量是可变的，故模型二中的重物球为迭代方程的新增自变量。根据（5.10）、（5.11）可以得到：

???

???????????-+--+---αθαπραθcos cos )2025.1(cos )sin 2(25.1cos min 11221621T T v y T v g T h T s 风风海（5.12）

??≤≤≤≤≤≤10000

12006000090

0..11q m T t s θ （5.13）其中增加了重物球质量的变量和对它的约束。在（5.12）中的量由（5.14）

求得

)

,,(],,,[11216216216216q y x m T f T T y x θ= （5.14）

5.2.2模型二的求解

问题二中第一问的假设是在第一问的基础上得到的故可以使用第一问的模型进行求解，通过MATLAB 实现得到:

风速36m/s 时26.33391=θ，浮标吃水0.7482米，游动区域为以海面上锚纵坐标为圆心半径为18.9578米的圆形区域，从下往上四节钢管在竖直方向倾斜角度为8.4463度、8.4225度、8.3989度、8.3753度。钢桶在竖直方向倾斜角度

为

8.5294=桶θ。锚链形状如图10所示

图10

由此可见在风速36km/h 的情况下，重物球的质量为1200kg 的时候，

26.3339

1=θ，8.5294=桶θ。根据题目二中161≤θ和5≤桶θ的要求，可以得知重物球的当前质量是不能满足要求的。所以需要求解（5.12）、（5.13）。

在（5.12）中q m 按10的步长在]10000,1200

[范围依次取不同的q m 值，每个对应的取值都有对应的浮标吃水深度h 、锚链末端与锚的链接处与海床的夹角

1θ、钢桶的倾斜角度桶θ。截取]6000,2000[部分的图像如图11所示

图11

通过条件：

??≤≤≤≤<<5

016021桶θθh h s （5.15）来筛选出对应可行的q m 值。通过MATLAB 软件实现得到重物球的质量满足条

件的一组为3700到5320中间按10个间隔，通过图像可以看出各个量的值在]5320,3700[无巨大波动。考虑三个限制条件的数据变化缓慢平滑，说明在]5320,3700[范围内可以取极小步长或连续取值，故可以得到重物球的取值范围为最低3700千克，最大5320千克，可以使得钢桶的倾斜角度不超过5度，锚链在锚点与海床的夹角不超过16度。

5.3问题三

5.3.1模型三的建立

根据问题三的叙述，和问题分析中的分析可以得知水深、风速、水速为变量。题目中给出水深的变化范围，水速和风速的变化范围，在上述范围内不同取值的组合构成了不同的情况；和模型二相比锚链的型号、锚链的长度、重物球的质量都变成了变量。需要优化的变量为吃水深度、游标游动区域和钢桶的倾斜角。对本题题意的理解可以解释为：设计一种系泊系统也就是确定锚链型号、锚链长度、重物球的量使得在不同水深、风速、水速情况下满足要求且最优。由此可以建立优化模型：

???

???????∑∑∑t t t x h max min 桶θ （5.16）

?????

=====≥<≤≤≤><<≤≤06

.5,925.2,5.1,735.0,249.018.0,15.0,12.0,105.0,078.020,18,1636,12,05

.1,5.0,00

10000016002050..11max l q m l h v v m T x h t s 海风水桶θθ （5.17） ),,,,,,,,(],,,[11216216216216q l y x m n m l h v v T f T T y x 海风水θ= （5.18）

根据水深、海水速度、风速的取值范围取得s 种不同的情况：

??????=20365.118125.01600

s 矩阵s 的1,2,3列为风速、水流速、水深的3个取值。每一行可以当做一种

情况记为t 。

对于（5.18）式所涉及的迭代模型，由于引入了水流速的作用所以产生水流力，再次进行受力分析如图12：

图12

钢桶的为：

g m T T t yn n y 2

)1(+=+

n xn n x v l T T θsin 1.562

)1(水桶-=+ （5.19）

钢管的为：

g m T T g yn n y 2

)1(+

=+ n xn n x v l T T θsin 35.92

)1(水管-=+ （5.20）

其余各式的同（5.2）、（5.3）式。

5.3.2模型三的求解

通过MATLAB 编程求解得到最优设计：选用Ⅲ型锚链，锚链长度为27.24m ，

重物球质量为2580kg 。选用此设计方案，如果水速、风速为0，水深为16m 时得到：吃水0.6m ，游动区域0m ，各节为90度，钢桶倾角为0度。

如果水速、风速最大，水深为18m 时，得到：吃水0.82m ，游动区域23.2m ，钢桶倾角为4.47度。

六、模型检验

6.1模型求解算法验证

对于问题一二三中通过迭代关系确立的表达式的优化求解是相同的，在所有模型中使用的都是多层搜索算法来寻优。由于该算法对初始精度取值敏感，一开始就进入局部最优搜索，所以对于问题一的求解，不能确定其精确性。这里对模型一的结果使用遗传算法进行检验，遗传算法是一种基于概率的进化算法，其特点是对结构对象操作，不用求导和函数连续性的限定，有隐并行性和较好的全局寻优能力；能自动获取优化的搜索空间，自适应能力可以自发调整搜索方向，不用给定的规则。

一般遗传算法的步骤如下：

图13

这里通过MATLAB软件中的ga函数实现遗传算法（代码见附录），由于遗传算法是基于概率的算法所以，每次的结果都会不同，这里选取多组结果如下表：锚末端与海床夹角浮标本身吃水深度(m) 浮动范围(m) 钢桶与竖直平面夹角

10.4001 0.7193 16.5369 2.2895

1.01 0.7166 16.4678 1.872

0 0.7164 16.2684 1.7802

10.4001 0.7193 16.5369 2.2895

1.01 0.7166 16.4678 1.872

0 0.7164 15.2868 1.5219

4.5715 0.7175 16.522 2.0184

0.7564 0.7166 16.4627 1.8621

13.1314 0.7204 16.5091 2.4318

0 0.7164 15.5756 1.5928

0.7182 0.7166 16.4619 1.8606

0 0.7164 15.9514 1.6912

10.3339 0.7193 16.5373 2.2861

9.2483 0.7189 16.5417 2.2325

0.0001 0.7164 16.4467 1.8328

13.3973 0.7205 16.5051 2.4462

4.1864 0.7174 16.5177 2.0019

0.006 0.7164 16.4468 1.8331

0.4533 0.7165 16.4564 1.8503

9.7874 0.7191 16.54 2.2589

4.47052 0.71765 16.35195 1.99119

最后一行为各个量的平均值，通过和多层搜索算法结果的对比发现，两种算法得到的结果相差不大，考虑到遗传算法的随机性，在真正利用算法求解模型时结果的波动太大，得到的结果往往不是很好；初始高精度的多层搜索算法，虽然运行缓慢，代码效率不高，但是结果固定，在准确度上也较高。

6.2模型仿真验证

对于模型本身的验证，对于本题目涉及到的问题可以选用计算机仿真，也可以选用实物仿真。通过调查发现，对于复杂的系泊系统缺少相关的仿真软件，虽可以自己编写仿真程序，但是由于时间限制，自己编写仿真的程序在短时间内也无法完成，故这里在仿真验证上采用实物仿真的方式。

实物仿真如下图14所示，观测数据的细节见附录，各个部分的数据分别为

图14

锚的质量0.16kg ，浮标直径0.06m ，浮标高0.06m ，浮标质量0.017kg ，钢管长度0.029m ，钢桶长度0.029m ，钢桶质量0.02kg ，钢桶体积3510653.1m -?，重物球质量0.02kg ，锚链长度0.09m ，锚链节数8节，锚链质量0.01kg ，水深0.17m ，吃水0.015m 。

根据以上所给的数据带入模型一可以得到一系列结果，考虑到测量的难度和准确度，这里取1θ和桶θ作为验证模型的参数。从附件中的细节测量图可以测得，实际的 9.91=θ， 3.40=桶θ。然而模型求得的 2.131=θ， 7.33=桶θ。可以看出实物仿真的结果和模型求解结果误差，误差的产生主要由两方面的原因：

1.实物的基本参数测量的误差，以及仿真结果测量的误差

2.模型对于现实问题存在简化，模型存在不精确的地方

七、模型的评价与改进

7.1模型的评价

在实际中锚链是分节的，该模型考虑到了实际情况下锚链的离散特性，没有把锚链部分当做连续的垂链线处理，使得该模型在精确度上比把锚链当做连续分析要好。

对锚链的每节进行分析，然后倒推得到的一个迭代表达式，然后优化找到合适的值。这样该模型的一些参数是靠优化求得。如果该模型有解析解显然要优于优化结果。

7.2模型的改进

对于浮标系统部分，在模型一二三中受力分析简化为：

实际的受力分析为：

改进后的模型更为贴近实际，但是在求解时难度会增加。

八、参考文献

[1] 郑瑞杰. 锚泊系统受力分析[D]. 大连: 大连理工大学, 2006.

[2] 王丹. 一般状态下悬链线方程的应用[J]. 武汉: 武汉理工大学交通学院, 2007: 26-28

[3] 乔东升,欧进萍. 深水悬链锚泊系统静力分析[J]. 航船工程, 2009, 38(2): 120-124.

[4] 潘甜. 组合锚泊系统的分析与研究[D]. 武汉: 武汉理工大学, 2010.

[5] 王艳妮. 海洋工程锚泊系统的分析研究[D]. 哈尔滨:哈尔滨工程大学, 2006.

[6] 肖越, 王言英. 浮体锚泊系统计算分析[J]. 船舶工程, 2005, 45(5): 682-686.

[7] 陈云鹤.赵晶.考虑水流力作用的锚链力计算方法[J].野战工程专

业,2016,17(2):116-120.

附录

问题一主程序q1_201666.m

注：此段为风速为12m/s的程序，精度为 p_jingdu=100，24改相应数据即可，即：v=24, p_jingdu=130

clc;clear

v=12;

rou=1025;

vfu=pi*2;

hself=1000/(rou*pi);

hfeng=2-hself;

%F=0.625*s*v^2;

ll=1;

Vt=0.15^2*ll*pi;

Vg=0.025^2*ll*pi;

mt=100-Vt*rou;

mg=10-Vg*rou;

mq=1200;

l=0.105;

J=[0,0];

ml=7*l;

g=9.8;

%初始变量范围

h_size=[0,8000];

w_size=[0,30];

%开始进行多层搜索

cengshu=5;%搜索层数

for cishu=1:cengshu

p_jingdu=100;

%每层精度处理

jingdu_h=8000/((p_jingdu)*((cishu-1)^3+1));

jingdu_w=30/((p_jingdu)*((cishu-1)^3+1));

%搜索循环几个变量几层循环

i=1;

for T0=h_size(1):jingdu_h:h_size(2)

for theta=w_size(1):jingdu_w:w_size(2)

e2(i)=q1_2016666_duoceng([T0,theta,mq,v]); note(i,:)=[T0,theta];

i=i+1;

end

emin=min(abs(e2));%m每层e的最小值

good=note(find(abs(e2)==emin),:);%最优解

%设置下层搜索范围为本层搜索的最优解值+-本层精度值 h_size=[good(1)-jingdu_h,good(1)+jingdu_h];

w_size=[good(2)-jingdu_w,good(2)+jingdu_w];

end

%plot(abs(e2),'.')

thetabest=good(1,2)

Tbest=good(1,1);

T(1,1)=cosd(thetabest)*Tbest;

T(1,2)=sind(thetabest)*Tbest;

theta0(1)=thetabest;

suml=22.05/l;

for n=1:suml-1

T(n+1,1)=T(n,1);

T(n+1,2)=T(n,2)+0.5*ml*g;

theta0(n)=atan((T(n,2)/T(1,1)))*180/pi;

J(n+1,1)=J(n,1)+l*cosd(theta0(n));

J(n+1,2)=J(n,2)+l*sind(theta0(n));

end

%211个节点

n=n+1;

T(n+1,1)=T(n,1);

T(n+1,2)=T(n,2)+0.5*ml*g+mq*g;

theta0(n)=atan((T(n,2)/T(1,1)))*180/pi;

J(n+1,1)=J(n,1)+l*cosd(theta0(n));

J(n+1,2)=J(n,2)+l*sind(theta0(n));

%钢桶

n=n+1;

T(n+1,1)=T(n,1);

T(n+1,2)=T(n,2)+0.5*mt*g;

theta0(n)=atan((T(n,2)/T(1,1)))*180/pi;

J(n+1,1)=J(n,1)+ll*cosd(theta0(n));

J(n+1,2)=J(n,2)+ll*sind(theta0(n));

%钢管

for m=n+1:n+4

T(m+1,1)=T(m,1);

T(m+1,2)=T(m,2)+0.5*mg*g;

theta0(m)=atan(T(m,2)/T(1,1))*180/pi;

J(m+1,1)=J(m,1)+ll*cosd(theta0(m));

J(m+1,2)=J(m,2)+ll*sind(theta0(m));

end

h=18-J(end,2)

h1=hself+sqrt(T(end,1)^2+T(end,2)^2)*sind(theta0(end))/(g*pi*rou)

axismove=J(end,1)

gejie=theta0(end-3:end)

tongjiao=90-theta0(end-4)

plot(J(1:end-5,1),J(1:end-5,2))

xlabel('长度/米');ylabel('长度/米');title('风速12m/s时锚链的形状')

问题一函数q1_2016666_duoceng.m

function e2=q1_2016_duoceng(p)

g=9.8;

rou=1025;

hself=1000/(rou*pi);

ll=1;

Vt=0.15^2*ll*pi;

Vg=0.025^2*ll*pi;

mt=100-Vt*rou;

mg=10-Vg*rou;

l=0.105;

J=[0,0];

ml=7*l;

T(1,1)=cosd(p(2))*p(1);

T(1,2)=sind(p(2))*p(1);

theta0(1)=p(2);

suml=22.05/l;

for n=1:suml-1

T(n+1,1)=T(n,1);

T(n+1,2)=T(n,2)+0.5*ml*g;

theta0(n)=atan((T(n,2)/T(1,1)))*180/pi;

J(n+1,1)=J(n,1)+l*cosd(theta0(n));

J(n+1,2)=J(n,2)+l*sind(theta0(n));

end

%211个节点

n=n+1;

T(n+1,1)=T(n,1);

T(n+1,2)=T(n,2)+0.5*ml*g+p(3)*g;

theta0(n)=atan((T(n,2)/T(1,1)))*180/pi;

J(n+1,1)=J(n,1)+l*cosd(theta0(n));

J(n+1,2)=J(n,2)+l*sind(theta0(n));

%钢桶

n=n+1;

T(n+1,1)=T(n,1);

T(n+1,2)=T(n,2)+0.5*mt*g;

theta0(n)=atan((T(n,2)/T(1,1)))*180/pi;

J(n+1,1)=J(n,1)+ll*cosd(theta0(n));

J(n+1,2)=J(n,2)+ll*sind(theta0(n));

%钢管

for m=n+1:n+4

T(m+1,1)=T(m,1);

T(m+1,2)=T(m,2)+0.5*mg*g;

theta0(m)=atan(T(m,2)/T(1,1))*180/pi;

J(m+1,1)=J(m,1)+ll*cosd(theta0(m));

J(m+1,2)=J(m,2)+ll*sind(theta0(m));

end

test1=T(end,1);

test2=1.25*(2-hself-T(end,2)/(pi*rou*g))*p(4)^2;

test3=(1.25*J(end,2)-20)*p(4)^2;

e2=abs(test1-test2)+abs(test1-test3);

问题二主程序q2_2016.m

注：第一小问风速为36m/s，修改问题一主程序即可，即：v=36, p_jingdu=130；第二小问主程序如下：

clc;clear

v=36;

rou=1025;

vfu=pi*2;

hself=1000/(rou*pi);

hfeng=2-hself;

ll=1;

Vt=0.15^2*ll*pi;

Vg=0.025^2*ll*pi;

mt=100-Vt*rou;

mg=10-Vg*rou;

mq=1200;

l=0.105;

J=[0,0];

ml=7*l;

g=9.8;

%初始变量范围

h_size=[0,6000];

w_size=[0,16];

o=1;

for mq=1200:10:7000

%开始进行多层搜索

cengshu=5;%搜索层数

for cishu=1:cengshu

p_jingdu=130;

%每层精度处理

jingdu_h=8000/((p_jingdu)*((cishu-1)^3+1));

jingdu_w=16/((p_jingdu)*((cishu-1)^3+1));

%搜索循环几个变量几层循环

i=1;

for T0=h_size(1):jingdu_h:h_size(2)

for theta=w_size(1):jingdu_w:w_size(2)

e2(i)=q1_2016666_duoceng([T0,theta,mq,v]); note(i,:)=[T0,theta];

i=i+1;

end

emin=min(abs(e2));%m每层e的最小值

good=note(find(abs(e2)==emin),:);%最优解

%设置下层搜索范围为本层搜索的最优解值+-本层精度值 h_size=[good(1)-jingdu_h,good(1)+jingdu_h];

w_size=[good(2)-jingdu_w,good(2)+jingdu_w];

end

thetabest(o)=good(1,2);

Tbest(o)=good(1,1);

T(1,1)=cosd(thetabest(o))*Tbest(o);

T(1,2)=sind(thetabest(o))*Tbest(o);

theta0(1)=thetabest(o);

suml=22.05/l;

for n=1:suml-1