2016数学建模国赛A题程序（原创）作者cclplus

杨利霞

2016数学建模国赛A题程序（原创）作者cclplus

2 H1 a+ H2 r; X1 c
9 r0 @3 f0 q* _
clear all;
close all;
clc
: L. X# j. E& o; yformat long
syms h S Fw Ff Ff1 a b c d l L F depth n pl m x1 y1 y t distance n a1 b1;
F=[];
theta=[];
v=24; %风速
, B& J/ M& ~* Wl=105*10^(-3); %锚链每节链环的长度
L=22.05; %锚链的总长度
num=0; %通过更改不在海床上的链节数得到一个最优解
" O* q) h; m, l: M. k' Anum1=round(L/l);
( g! y0 B( ]3 [num2=0;
lin=0/180*pi; %第一个链节与水平方向的夹角
4 _$ J2 a( n7 x/ D6 F7 {lin1=90/180*pi;
! [. g; L4 j1 u: k/ i' t- J9 @, v$ Xlin2=0;
m2=1200; %重物球质量
$ W& B0 g- \0 ?! Dpg=7.7*10^3; %重物球的密度（单位：kg/m^3)
depth=20; %水深
* O* I& O( o6 }' j% o! mpl=7; %锚链单位长度的质量
  Z0 w+ K6 L* d8 f  Wvh=0; %海水流速
1 y( q% k8 `% @$ tg=9.8; %可通过改变此语句来修改重力加速度，单位为m/s^2
p=1.025*10^3; %海水密度
2 b* W: T$ ]* T; f: h& P) [( N. W% ]M=1000; %浮标质量
4 @0 `' }* t/ X. l0 u2 Rm=10; %钢管质量
m1=100; %设备和钢桶总质量
y=0;
6 F$ x8 d( S4 [$ s7 T6 C; X+ a) Fd=1;
j1=0;
j2=0;
while(abs(y-d)>0.005)%在这里选择所需要的精度，
if (y>d)&&(num<round(L/l))
. G  }$ |* Z! J. O- Unum1=num;
num=round((num1+num2)/2);
+ g1 n  x2 W3 a/ J- _4 b% Oelseif (y<d)&&(num<round(L/l));
num2=num;
num=round((num1+num2)/2);
elseif (y<d)&&(num==round(L/l))
lin2=lin;
lin=(lin1+lin2)/2;
elseif(y>d)&&(num==round(L/l))
- h& {1 A" T: y: R2 u# Glin1=lin;
  x' H) Y0 s! z, l  i) E9 T) xlin=(lin1+lin2)/2;
/ @5 a6 M# y9 q9 B/ W  Y" ^end
/ f: ]) X$ |/ v, Z, K2 J7 d# h%钢桶受到的浮力
Ff1=p*g*pi*(0.3/2)^2;
  x4 d: \8 V  J6 C. l: u%钢管收到的浮力
Ff2=p*g*pi*(0.05/2)^2;
7 F6 j1 _/ K$ |/ g4 h%重物球所受浮力
Ffg=p*g*m2/pg;
. n% _, S% v" L# ^1 M6 e5 w$ r! t%重物球所受海水水流力
& R6 e+ o) C" u0 a& q5 L* iFhg=374*pi*((m2/pg/3/4)^(1/3))^2*vh^2;
%风对浮标受力面的投影面积
S=2*(2-h);
3 [+ _2 T* _+ c: a5 ?: j+ g: E( T%风对浮标产生的力
9 t0 d# I8 D# k; d$ k- f; RFw=0.625*S*v^2;
- L  _2 s" Z# H, S) S" S%浮标在水中的体积
V=pi*(2/2)^2*h;
%浮标所受到的浮力
6 {; W. |7 O# S5 X1 ?Ff=p*g*V;
; I( A7 U$ f1 O1 \- f* o5 A0 `( `%浮标受到海水的近似水流力
& R( O; c9 [0 r$ U. JFb=374*2*h*vh^2;
%钢桶受到海水的近似水流力
Fs1=374*0.3*vh^2;
, _, b* S8 g  M%钢管受到海水的水流力的近似值
( ]" ?# u% x! @' QFs=374*0.05*vh^2;
+ k, L! \/ P5 j9 M. B3 r# l%浮标浸没水中的高度
if num==round(L/l)
" ?# k/ M( b, H0 z4 uh=(m2*g+M*g+4*m*g+m1*g-Ff1-4*Ff2-Ffg+pl*L*g+(Fhg+4*Fs+Fs1)*tan(lin))/(p*g*pi-(1.25*v^2+374*vh^2)*tan(lin));
5 g1 c# [# ~! v( felse
( `9 w! {* @/ K8 J6 p' u3 `h=(m2*g+M*g+4*m*g+m1*g-Ff1-4*Ff2-Ffg+num*pl*l*g)/(p*g*pi);
end
/ l0 R0 L7 X' a0 T5 ^. M, Wa=Fw+Fb;
+ m( I8 J) m$ f2 L5 Mb=-M*g+Ff+(Fw+Fb)*tan(lin);
. _. m- o7 q& Nif j1==0
a=eval(a);
; O/ c! w+ J3 Q6 X  hb=eval(b);
else
end
4 N& c1 s$ N& I5 Y) }2 @F(1)=sqrt(a^2+b^2);
theta(1)=atan(b/a);
n=0;
+ k0 `, Y, e8 G- Zfor i=1:4
6 B9 `' v7 o$ Z%钢管受到海水的水流力
Fh(i)=374*0.05*sin(theta(i));
' i  K% B+ }1 j7 ^$ a+ Jn=n+Fh(i);
a=Fw+Fb+n;
# k4 R, ~% A4 k, Y' }- A" ?1 Kif j1==0
a=eval(a);
/ o" @9 j2 x+ @6 c' f! k& xelse
end
' `' t" x6 |( F. c0 Ab=F(i)*sin(theta(i))+p*g*pi*(50*10^(-3)/2)^2-m*g;
F(i+1)=sqrt(a^2+b^2);
* |: f$ V* V% w% c8 K7 r6 Ctheta(i+1)=atan(b/a);
3 Z9 r# x/ d" K. f  Zend
4 a! ^, K, h" P( Z7 ^7 \c=0;
& M9 b3 J  r# l) G5 Z  \for i=1:5
c=c+sin(theta(i));
6 \2 u" P+ J, X1 Q* u) uend
d=depth-c-h;
y1=lin;
( G( l" z0 a* z2 ^4 i! Gdistance=0;
# k/ ~! M) `4 [% lif num==round(L/l)
y=l*sin(y1);
* {  o/ \: s, h8 ^# ?+ }$ }x1=Fw/sqrt(1-(sin(y1))^2);
+ m3 N( N1 n2 L+ v; y0 c2 f% `for i=1:num-1
" g. D7 B" v/ R3 z1 q! y% nm=(x1*sin(y1)+i*pl*l*g)/sqrt((x1*sin(y1)+i*pl*l*g)^2+Fw^2);
m=m*l;
" E" Y7 }/ B4 j- O2 n' b7 C. oy=y+m;
+ n) @! k( A$ ~$ S5 Bn=Fw/sqrt((x1*y1+i*pl*l*g)^2+Fw^2)*l;
if j1==0
. E" ?( W- e* o) |9 s  w2 gn=eval(n);
else
end
& m2 S% }: q) f' A) N* n. ]0 Bdistance=distance+n;
if j1==0
y=eval(y);
else
- \' U; T9 {; F5 z$ r9 f% Tend
! i3 f3 u* L- |3 L, }  Oend
* M! u, Q" f  ~# n& I- ?( C& m8 eelse
y=y1*l;
distance=(round(L/l)-num)*l;
for i=1:num
x1=Fw/sqrt(1-(sin(y1))^2);
: }% Q& N. S2 l6 l& bm=(x1*y1+i*pl*l*g)/sqrt((x1*y1+i*pl*l*g)^2+Fw^2)*l;
$ K; I) z& `; ?2 py=y+m;
n=Fw/sqrt((x1*y1+i*pl*l*g)^2+Fw^2)*l;
' f/ [$ v3 t  Q; S# d2 @if j1==0
y=eval(y);
n=eval(n);
else
end
7 ^3 a( A* X& C# q( adistance=distance+n;
end
# @3 z  y  _; Z7 W9 ~1 pend
1 o4 I$ ]( L' Fm=0;
j1=1;
7 V, f! Z' x+ _: o, U& A! \j2=j2+1;
end
%钢桶受到的浮力
* J7 W' M& M! H# p- GFf1=p*g*pi*(0.3/2)^2;
; _8 P& s( X5 m  }! b%钢管收到的浮力
Ff2=p*g*pi*(0.05/2)^2;
, `) I1 ]) T: S+ I%重物球所受浮力
Ffg=p*g*m2/pg;
/ `$ b, j% j* q4 U; v%重物球所受海水水流力
Fhg=374*pi*((m2/pg/3/4)^(1/3))^2*vh^2;
%风对浮标受力面的投影面积
S=2*(2-h);
%风对浮标产生的力
9 s: D* C9 {3 @5 j1 p. e9 j0 O! W3 EFw=0.625*S*v^2;
% N+ j% d" p: R, @- W$ i' m  V%浮标在水中的体积
  ?* t* N5 V( |4 L, YV=pi*(2/2)^2*h;
. q3 |+ S% }* A1 Q( p%浮标所受到的浮力
! G1 p6 X2 N+ b/ f) ?5 J! C2 y. DFf=p*g*V;
%浮标受到海水的近似水流力
8 d0 S4 E9 {* k6 T' rFb=374*2*h*vh^2;
%钢桶受到海水的近似水流力
) a" w+ i# s: {7 x1 X2 S3 pFs1=374*0.3*vh^2;
%钢管受到海水的水流力的近似值
4 T& X) ?5 V/ _' w5 L$ l) v, qFs=374*0.05*vh^2;
%浮标浸没水中的高度
if num==round(L/l)
h=(m2*g+M*g+4*m*g+m1*g-Ff1-4*Ff2-Ffg+pl*L*g+(Fhg+4*Fs+Fs1)*tan(lin))/(p*g*pi-(1.25*v^2+374*vh^2)*tan(lin));
else
$ q3 |, N3 I) k- {% m, mh=(m2*g+M*g+4*m*g+m1*g-Ff1-4*Ff2-Ffg+num*pl*l*g)/(p*g*pi);
6 n! C# B' a; m' p$ h# F' S3 Gend
a=Fw+Fb;
b=-M*g+Ff+(Fw+Fb)*tan(lin);
F(1)=sqrt(a^2+b^2);
theta(1)=atan(b/a);
. N% B) q( ?; `* @6 sn=0;
for i=1:4
%钢管受到海水的水流力
8 h6 T3 c7 p: I, _! ]) W5 pFh(i)=374*0.05*sin(theta(i));
n=n+Fh(i);
" `/ g4 d# Y  e8 e: s+ a% @  Qa=Fw+Fb+n;
2 {0 P1 w1 F% V/ X( y2 |" a& ?b=F(i)*sin(theta(i))+p*g*pi*(50*10^(-3)/2)^2-m*g;
F(i+1)=sqrt(a^2+b^2);
, Y' f2 N: {9 |$ f7 {2 ctheta(i+1)=atan(b/a);
' i& r7 n- L7 A" Y5 ]- ^1 ^end
, y" B- G0 Z) j& S9 Hdisp('输出钢管和钢桶的倾斜角度（角度制）')
8 Z3 r8 `/ ]0 lth=90-theta*180/pi
m=85*pi/180;
+ B% a, C  z: H, q, B. t* D2 C% J: gif theta(5)>m
disp('钢桶的倾斜角足够小，测量准确')
1 s. y; U% y& e# welse
) z* P% q4 M9 a2 X- ]disp('钢桶的倾斜角过大')
  {6 E! B9 g; P; }# iend
c=0;
6 l& |3 v) C1 U4 Xfor i=1:5
% t" a8 }" c0 G, q9 cc=c+sin(theta(i));
- y; n2 r) p6 Q- m+ ?! M. i: {end
d=depth-c-h;
7 X+ F( y; ?7 k* {! d1 Q- [y1=lin;
" x3 l8 N6 F, d$ Kdistance=0;
if num==round(L/l)
. |7 P: q4 ?' O+ ry=l*sin(y1);
' n0 R- S7 ?4 v& C5 y( G( Qx1=Fw/sqrt(1-(sin(y1))^2);
for i=1:num-1
( F* \7 D$ k/ X3 J% q+ x% fm=(x1*sin(y1)+i*pl*l*g)/sqrt((x1*sin(y1)+i*pl*l*g)^2+Fw^2);
m=m*l;
y=y+m;
n=Fw/sqrt((x1*y1+i*pl*l*g)^2+Fw^2)*l;
distance=distance+n;
plot(distance,y,'o')
$ Q# t" t' m) Q' t% b( i" W6 xhold on
5 ^6 R* C8 Q/ q/ @" ~end
else
# u* I0 \# \6 {, My=y1*l;
9 a+ C' W2 _# K: O# k5 \& dfor i=1:round(L/l)-num
distance=i*l;
y=0;
5 A8 |# ?9 p9 Z) @plot(distance,y,'o')
. C$ }1 C* o/ p* t1 C& Dhold on
grid on
9 M, d, _' Q" d9 t) t( o1 L1 d' W5 Qend
" Z0 \# W. z8 n) nfor i=1:num
+ W$ y) {# g( r. cx1=Fw/sqrt(1-(sin(y1))^2);
8 S/ G& r5 e" Ym=(x1*y1+i*pl*l*g)/sqrt((x1*y1+i*pl*l*g)^2+Fw^2)*l;
  T4 u+ D5 [3 K# e7 cy=y+m;
8 U# q; K/ y: a3 {0 }n=Fw/sqrt((x1*y1+i*pl*l*g)^2+Fw^2)*l;
if j1==0
y=eval(y);
n=eval(n);
8 I0 Y! [: F) \9 I2 \else
8 w8 D' h) \+ m# Q( U! u8 P) B  uend
) g% A  L8 @7 v! B$ S5 M; Kdistance=distance+n;
. o+ }  P: z5 u1 Qplot(distance,y,'o')
hold on
end
# ?; q2 V4 `3 p& jend
0 b: M# O, T& g1 Pm=0;
for i=1:5
m=m+cos(theta(i));
end
0 v& L& g8 [/ w  S& H: \0 d5 P%浮标的运动半径
disp('输出浮标的运动半径')
: ~- O" e" U" `9 Fans=distance+m
, y' q6 d1 b! M* c


) l/ e7 o1 L% y( b7 V5 e
/ H4 y: }& _9 f/ y


: p; |! [, _& P

571334077

2333333333333333333333