数学建模初步

一只小船渡过宽为d的河流，目标是起点A正对着的另一岸B点，已知河水流速v1与船在静水中的速度v2之比为k，

1）建立小船航线的方程，求其解析解；

2）设d=100,v1=1m/s,v2=2m/s，用数值解法求渡河所需时间、任意时刻小船的位置及航行曲线，作图，并与解析解比较。

3）若流速v1为0，0.5，1.5，2 (m/s)，结果将如何.

小船渡河问题

一、模型准备

首先，我们建模的目的就是确定小船在不同的水速以及小船在静水中的不同船速下渡河的航线问题。其次，我们要明确影响小船航线变化的因素有哪些?最明显的因素是河水流速、小船在静水中的船速以及小船行驶的方向。最后，我们就可以建立小船航线与河水流速、小船在静水中的船速以及小船行驶的方向之间的数学关系式。

二、模型假设

1．船行驶路线为连续曲线，终点为起点A对面的点B。

2．船在行驶过程中始终向着B点前进，即船速

v始终指向点B；并且过程中船

2

的合速度v始终与船的行驶轨迹相切。

3．该段河流为理想的直段，水速

v与河岸平行。

1

三、模型建立

该问题模型为典型的微分方程模型。问题中船在静水中的速度及流水速度是不变的，且有简单的比例关系，容易得出船的航线解析方程。同时船的航线的起点和终点已经确定，这为微分方程求数值解提供了初值条件。为了刻画小船在不同时刻的位置，我们必须先建立直角坐标系。我们以船的终点B点为坐标原点，水速

v的方向为x轴正方向，以B到A的方向为y轴正方向，建立直角坐标系。

1

如下图：

B x轴

v

2

v

1

A

y轴

为了方便起见，我们把所用到得符号意义归纳一下，如下： d 河宽； 1v 河水流速；

2v 小船在静水中的速度；

k 河水流速1v 与船在静水中的速度2v 之比为k ； θ 小船在静水中的速度2v 的反向与y 轴正向的夹角； )(t r 在t 时刻时，小船到坐标原点B 的距离； )(t x 在t 时刻时，船在x 轴方向上的位移为x ； )(t y 在t 时刻时，船在y 轴方向上的位移为y ；

x v 或

dt dx

在t 时刻，船在x 轴方向上的合速度； y v 或dt

dy

在t 时刻，船在y 轴方向上的合速度；

在t 时刻，将船在静水中的速度2v 分别在x 轴、y 轴方向上的进行分解，可得： θsin *21v v v x -= ----------(1)式 θcos *2v v y -= ----------(2)式 由假设(2),船在静水中的速度2v 始终指向点B ,所以可得： y

x =

θtan 进而可求出 2

2

sin y

x x +=θ，2

2

cos y

x y +=θ ----------(3)式

将(3)式分别带入(1)、(2)式，得：

2221*y

x x v v v dt dx x +-==

222*y

x y

v v dt dy y +-==

上述两式相除，可得：

x y x k y

y

x x k y x y

y x x v v y x y

v dx dy -+-=+-+-

=+-+-=222

2222221222***

————(4)式

再对(4)式两边同时积分即可求出小船的航线，然而实际用Matlab 计算时，求不

出其解析解，而问题让我们求其解析解，故我们不是直接对(4)式积分，转换下形式。

我们令θsin *r x =，θcos *r y =将直角坐标化为极坐标，由导数的的链式法则，我们可以得到：

θθθd r dr dx *cos **sin += ————(5)式 θθθd r dr dy *sin **cos -= ————(6)式 因为dt dx v x =，dt

dy

v y =，将其带入(5)、(6)两式，可得：

θθθd r dr dt v x *cos **sin *+= θθθd r dr dt v y *sin **cos *-= 将上述两式相除，再把(1)、(2)式带入得：

k

k r d dr --+=

)

tan *tan sin (cos *θθθθθ 两边积分可得r 关于θ的函数，再由θsin *r x =，θcos *r y =，进而可得出x 关于θ的函数、y 关于θ的函数，即的小船运动轨迹的参数方程。

四、模型求解

1、第一问

利用matlab 软件求解微分方程，在命令窗口中输入：

dsolve('Dy=y*(cos(x)+sin(x)*tan(x)-k*tan(x))/(-k)','y(0)=d','x')

(其中y 相当于r ，x 相当于θ)，可求出r 关于θ的函数为：

2

11)

+ )(sin(*1) - )(sin( i

*)k 2 i

*pi exp(-*)k 2 1) )log(sin( exp(-*)k 2 1) - )log(sin( exp(*d 2

1

θθθθ+=r

进而得到小船运动轨迹的参数方程为：

θθθθθsin *1)

+ )(sin(*1) - )(sin( i

*)k 2 i

*pi exp(-*)k 2 1) )log(sin( exp(-*)k 2 1) - )log(sin( exp(*d 2

12

1

+=

x

θθθθθcos *1)

+ )(sin(*1) - )(sin( i

*)k 2 i

*pi exp(-*)k 2 1) )log(sin( exp(-*)k 2 1) - )log(sin( exp(*d 2

121

+=

y

2、第二问

当100=d ，11=v ，22=v 时 程序如下：

建立m-文件fun.m 如下:

function Xdot=fun(t,x) v1=1; v2=2;

Xdot=zeros(2,1); if (norm(x)>1e-5)

Xdot(1)=v1-v2*x(1)/sqrt(x(1)^2+x(2)^2);%x(1)为小船的横坐标、x(2)为纵坐Xdot(2)=-v2*x(2)/sqrt(x(1)^2+x(2)^2); 标 else

Xdot=[0,0]; End

取t 0=0，t f =75，输入命令：

t0=0; tf=75;

x0=[0,100];

[t,x]=ode23('fun',[t0 tf],[0 100]);

[t,x(:,1),x(:,2)] %在t 时刻，小船的横坐标x(1)、纵坐标x(2)

plot(t,x(:,1),'g',t,x(:,2),'r');%绿线为小船的横坐标x(1)、红线为 title('小船渡河航线'); 纵坐标x(2) xlabel('时间t');

ylabel('小船航线横、纵坐标值x 、y');

求解得：

时间t 横坐标x 纵坐标y

0 0 100.0000 0.0001 0.0001 99.9998 0.0005 0.0005 99.9990 0.0025 0.0025 99.9950 0.0125 0.0125 99.9750 0.0625 0.0624 99.8750 0.3125 0.3115 99.3750 1.5625 1.5378 96.8752 6.0539 5.6723 87.9002 12.7098 10.9480 74.6636 20.2098 15.5068 59.9613 27.7098 18.3630 45.7078

35.2098 19.2574 32.2550 42.7098 17.9592 20.1285 47.5442 15.8835 13.3754 52.3786 12.8808 7.7266 55.2996 10.6502 4.9764

58.2206 8.1601 2.7804

59.9559 6.5766 1.7669

61.6913 4.9334 0.9745

62.7117 3.9445 0.6178

63.7321 2.9431 0.3406

64.3337 2.3481 0.2160

64.9352 1.7505 0.1193

65.2912 1.3958 0.0757 65.6473 1.0405 0.0419

65.8585 0.8295 0.0266

66.0698 0.6185 0.0147 66.1953 0.4930 0.0093 66.3207 0.3676 0.0052 66.3953 0.2931 0.0033 66.4699 0.2185 0.0018 66.5142 0.1742 0.0012 66.5585 0.1299 0.0006 66.5891 0.0993 0.0004 66.6197 0.0687 0.0002 66.6400 0.0484 0.0001 66.6603 0.0281 0.0000 66.6702 0.0181 0.0000 66.6802 0.0082 0.0000 66.6852 0.0031 0.0000 66.6878 0.0006 -0.0000 66.6881 0.0003 -0.0000 66.6882 0.0001 -0.0000 66.6883 0.0001 -0.0000 66.6884 0.0000 0.0000 66.6884 0.0000 0.0000 66.6884 0.0000 0.0000 66.6885 0.0000 0.0000 66.6892 0.0000 0.0000 66.6926 0.0000 0.0000 66.7093 0.0000 0.0000

66.7928 0.0000 0.0000

67.2107 0.0000 0.0000 69.3003 0.0000 0.0000 75.0000 0.0000 0.0000

当流速1v 为1时渡河航线如下图：

0102030

4050607080

-20

20

40

60

80

100

X: 66.71

Y: 1.724e-14

小船渡河航线

时间t

小船航线横、纵坐标值x 、y

从上图中可以看出，在流水速度为1时, 小船沿直线到达B 点的渡河时间为66.71秒。

3、第三问

当100=d ，12=v ，时

(1)若流速1v 为0时，我们用数值解法可以求出渡河所需时间、任意时刻小船的位置及航线的横、纵坐标。 程序如下：

建立m-文件fun.m 如下:

function Xdot=fun(t,x) v1=0; v2=2;

Xdot=zeros(2,1); if (norm(x)>1e-5)

Xdot(1)=v1-v2*x(1)/sqrt(x(1)^2+x(2)^2);%x(1)为小船的横坐标、x(2)为纵坐Xdot(2)=-v2*x(2)/sqrt(x(1)^2+x(2)^2); 标

else

Xdot=[0,0];

End

取t0=0，t f=70，输入命令：

t0=0;

tf=70;

x0=[0,100];

[t,x]=ode23('fun',[t0 tf],[0 100]);

[t,x(:,1),x(:,2)] %在t时刻，小船的横坐标x(1)、纵坐标x(2)

plot(t,x(:,1),'g',t,x(:,2),'r');%绿线为小船的横坐标x(1)、红线为title('小船渡河航线'); 纵坐标x(2)

xlabel('时间t');

ylabel('小船航线横、纵坐标值x、y');

求解得：

时间t横坐标x纵坐标y

0 0 100.0000

4.0000 0 92.0000

11.5000 0 77.0000

19.0000 0 62.0000

26.5000 0 47.0000

34.0000 0 32.0000

41.5000 0 17.0000

49.0000 0 2.0000

49.9375 0 0.1250

49.9961 0 0.0078

49.9998 0 0.0005

50.0000 0 0.0000

50.0000 0 0.0000

50.0000 0 0.0000

50.0000 0 0.0000

50.0000 0 0.0000

50.0002 0 0.0000

50.0011 0 0.0000

50.0053 0 0.0000

50.0266 0 0.0000

50.1328 0 0.0000

50.6640 0 0.0000

53.3199 0 0.0000 60.8199 0 0.0000 68.3199 0 0.0000 75.0000 0 0.0000

当流速1v 为0时渡河航线如下图：

10

20

30

4050

607080

01020304050607080

90100X: 50

Y: 6.676e-06

小船渡河航线

时间t

小船航线横、纵坐标值x 、y

从上图中可以看出，在静水中, 小船沿直线到达B 点(因为小船的横坐标始终为0), 渡河时间为50 秒,这与我们的直观经验是相符合的。

(2)若流速1v 为0.5时，我们用数值解法可以求出渡河所需时间、任意时刻小船的位置及航线的横、纵坐标。 程序如下：

建立m-文件fun.m 如下:

function Xdot=fun(t,x) v1=0.5; v2=2;

Xdot=zeros(2,1); if (norm(x)>1e-5)

Xdot(1)=v1-v2*x(1)/sqrt(x(1)^2+x(2)^2);%x(1)为小船的横坐标、x(2)为纵坐

Xdot(2)=-v2*x(2)/sqrt(x(1)^2+x(2)^2); 标

else

Xdot=[0,0];

end

取t0=0，t f=75，输入命令：

t0=0;

tf=75;

x0=[0,100];

[t,x]=ode23('fun',[t0 tf],[0 100]);

[t,x(:,1),x(:,2)] %在t时刻，小船的横坐标x(1)、纵坐标x(2)

plot(t,x(:,1),'g',t,x(:,2),'r');%绿线为小船的横坐标x(1)、红线为title('小船渡河航线'); 纵坐标x(2)

xlabel('时间t');

ylabel('小船航线横、纵坐标值x、y');

求解得：

时间t横坐标x纵坐标y

0 0 100.0000

0.0002 0.0001 99.9997

0.0010 0.0005 99.9981

0.0050 0.0025 99.9901

0.0250 0.0125 99.9501

0.1250 0.0624 99.7501

0.6250 0.3105 98.7501

3.1250 1.5126 93.7503

8.6298 3.9196 82.7465

15.6743 6.4600 68.6928

23.1536 8.3705 53.8486

30.2590 9.2543 39.8935

36.4877 9.0683 27.8845

41.6376 8.0156 18.2481

44.7168 6.8662 12.7094

47.1756 5.5859 8.4818

49.0117 4.3617 5.5004

50.3725 3.2658 3.4392

51.2482 2.4533 2.2134

51.8836 1.7989 1.3952

52.3108 1.3218 0.8925

52.6168 0.9581 0.5647

52.8277 0.6949 0.3599

52.9779 0.5001 0.2280

53.0828 0.3599 0.1450

53.1573 0.2578 0.0919 53.2098 0.1846 0.0583 53.2471 0.1319 0.0370 53.2734 0.0941 0.0235 53.2921 0.0671 0.0149 53.3053 0.0478 0.0094 53.3147 0.0341 0.0060 53.3214 0.0242 0.0038 53.3261 0.0172 0.0024 53.3295 0.0123 0.0015 53.3319 0.0087 0.0010 53.3336 0.0062 0.0006 53.3353 0.0036 0.0003 53.3361 0.0024 0.0002 53.3369 0.0011 0.0001 53.3373 0.0006 0.0000 53.3375 0.0003 0.0000 53.3376 0.0001 0.0000 53.3377 0.0000 -0.0000 53.3377 0.0000 -0.0000 53.3377 0.0000 -0.0000 53.3377 0.0000 -0.0000 53.3377 0.0000 -0.0000 53.3379 0.0000 -0.0000 53.3385 0.0000 -0.0000 53.3419 0.0000 -0.0000 53.3585 0.0000 -0.0000 53.4417 0.0000 -0.0000 53.8578 0.0000 -0.0000 55.9382 0.0000 -0.0000 63.4382 0.0000 -0.0000 70.9382 0.0000 -0.0000 75.0000 0.0000 -0.0000

当流速1v 为0.5时渡河航线如下图：

0102030

4050607080

-20

20

40

60

80

100

X: 53.34

Y: -6.093e-08

小船渡河航线

时间t

小船航线横、纵坐标值x 、y

从上图中可以看出，在流速1v 为0.5时, 小船的渡河时间为53.34秒。

(3)若流速1v 为1.5时，我们用数值解法可以求出渡河所需时间、任意时刻小船的位置及航线的横、纵坐标。 程序如下：

建立m-文件fun.m 如下:

function Xdot=fun(t,x) v1=1.5; v2=2;

Xdot=zeros(2,1); if (norm(x)>1e-5)

Xdot(1)=v1-v2*x(1)/sqrt(x(1)^2+x(2)^2);%x(1)为小船的横坐标、x(2)为纵坐Xdot(2)=-v2*x(2)/sqrt(x(1)^2+x(2)^2); 标 else

Xdot=[0,0]; End

取t0=0，t f=130，输入命令：

t0=0;

tf=130;

x0=[0,100];

[t,x]=ode23('fun',[t0 tf],[0 100]);

[t,x(:,1),x(:,2)] %在t时刻，小船的横坐标x(1)、纵坐标x(2)

plot(t,x(:,1),'g',t,x(:,2),'r');%绿线为小船的横坐标x(1)、红线为title('小船渡河航线'); 纵坐标x(2)

xlabel('时间t');

ylabel('小船航线横、纵坐标值x、y');

求解得：

时间t横坐标x纵坐标y

0 0 100.0000

0.0001 0.0001 99.9999

0.0003 0.0005 99.9994

0.0017 0.0025 99.9967

0.0083 0.0125 99.9834

0.0417 0.0625 99.9167

0.2083 0.3118 99.5834

1.0417 1.5461 97.9168

5.0076 7.1228 89.9950

11.4867 15.1009 77.1607

19.7419 23.0217 61.2825

29.7125 28.9684 43.5343

40.0237 31.0293 27.9248

49.1750 29.9683 17.1911

58.2611 27.1083 9.6035

67.6089 23.1249 4.6187

75.5215 19.3627 2.1732

80.8095 16.7637 1.2050

86.0976 14.1373 0.6003

89.4511 12.4649 0.3613

92.8047 10.7900 0.2017

95.4912 9.4474 0.1181

98.1777 8.1045 0.0635

99.6869 7.3499 0.0429

101.1961 6.5954 0.0278 103.2520 5.5674 0.0139 104.6651 4.8609 0.0081 106.0783 4.1543 0.0043 106.9615 3.7127 0.0027 107.8447 3.2711 0.0016 108.7312 2.8278 0.0009 109.6177 2.3846 0.0005 110.2264 2.0802 0.0003 110.8351 1.7759 0.0001 111.5960 1.3954 0.0000 112.3569 1.0150 0.0000 112.9926 0.6971 0.0000 113.6643 0.3613 -0.0000 113.9785 0.2042 0.0000 114.1356 0.1256 -0.0000 114.2927 0.0471 0.0000 114.3550 0.0159 -0.0000 114.3706 0.0081 0.0000 114.3784 0.0042 -0.0000 114.3823 0.0023 0.0000 114.3842 0.0013 -0.0000 114.3862 0.0003 0.0000 114.3864 0.0002 -0.0000 114.3867 0.0001 0.0000 114.3868 0.0000 -0.0000 114.3868 0.0000 0.0000 114.3868 0.0000 -0.0000 114.3869 0.0000 0.0000 114.3869 0.0000 0.0000 114.3870 0.0000 0.0000 114.3876 0.0000 0.0000 114.3907 0.0000 0.0000 114.4061 0.0000 0.0000 114.4833 0.0000 0.0000 114.8690 0.0000 0.0000 116.7977 0.0000 0.0000 126.4412 0.0000 0.0000 130.0000 0.0000 0.0000

当流速1v 为1.5时渡河航线如下图：

020406080100120140

-20

20

40

60

80

100

X: 114.3

Y: 1.206e-07

小船渡河航线

时间t

小船航线横、纵坐标值x 、y

从上图中可以看出，在流速1v 为1.5时, 小船的渡河时间为114.3秒。

(4)若流速1v 为2时，我们用数值解法可以求出渡河所需时间、任意时刻小船的位置及航线的横、纵坐标。 程序如下：

建立m-文件fun.m 如下:

function Xdot=fun(t,x) v1=2; v2=2;

Xdot=zeros(2,1); if (norm(x)>1e-5)

Xdot(1)=v1-v2*x(1)/sqrt(x(1)^2+x(2)^2);%x(1)为小船的横坐标、x(2)为纵坐Xdot(2)=-v2*x(2)/sqrt(x(1)^2+x(2)^2); 标 else

Xdot=[0,0]; End

取t0=0，t f=1100，输入命令：

t0=0;

tf=1100;

x0=[0,100];

[t,x]=ode23('fun',[t0 tf],[0 100]);

[t,x(:,1),x(:,2)] %在t时刻，小船的横坐标x(1)、纵坐标x(2)

plot(t,x(:,1),'g',t,x(:,2),'r');%绿线为小船的横坐标x(1)、红线为title('小船渡河航线'); 纵坐标x(2)

xlabel('时间t');

ylabel('小船航线横、纵坐标值x、y');

求解得：

时间t横坐标x纵坐标y

1.0e+03 *

0 0 0.1000

0.0000 0.0000 0.1000

0.0000 0.0000 0.1000

0.0000 0.0000 0.1000

0.0000 0.0000 0.1000

0.0000 0.0001 0.0999

0.0002 0.0003 0.0997

0.0008 0.0016 0.0984

0.0039 0.0075 0.0922

0.0101 0.0180 0.0800

0.0189 0.0299 0.0634

0.0300 0.0399 0.0449

0.0408 0.0453 0.0307

0.0508 0.0478 0.0211

0.0617 0.0491 0.0138

0.0721 0.0496 0.0091

0.0814 0.0498 0.0063

0.0902 0.0499 0.0044

0.0987 0.0500 0.0032

0.1071 0.0500 0.0023

0.1155 0.0500 0.0016

0.1238 0.0500 0.0012

0.1322 0.0500 0.0008

0.1405 0.0500 0.0006 0.1488 0.0500 0.0004 0.1572 0.0500 0.0003 0.1655 0.0500 0.0002 0.1738 0.0500 0.0002 0.1821 0.0500 0.0001 0.1904 0.0500 0.0001 0.1988 0.0500 0.0001 0.2071 0.0500 0.0000 0.2154 0.0500 0.0000 0.2237 0.0500 0.0000 0.2321 0.0500 0.0000 0.2404 0.0500 0.0000 0.2487 0.0500 0.0000 0.2570 0.0500 0.0000 0.2654 0.0500 0.0000 0.2737 0.0500 0.0000 0.2820 0.0500 0.0000 0.2903 0.0500 0.0000 0.2987 0.0500 0.0000 0.3070 0.0500 0.0000 0.3153 0.0500 0.0000 0.3244 0.0500 0.0000 0.3348 0.0500 0.0000 0.3469 0.0500 0.0000 0.3614 0.0500 0.0000 0.3796 0.0500 0.0000 0.4052 0.0500 0.0000 0.4428 0.0500 0.0000 0.4804 0.0500 0.0000 0.5798 0.0500 -0.0000 0.6555 0.0500 0.0000 0.7312 0.0500 -0.0000 0.8018 0.0500 0.0000 0.8593 0.0500 -0.0000 0.9144 0.0500 0.0000

0.9788 0.0500 -0.0000

1.0502 0.0500 0.0000 1.1000 0.0500 -0.0000

当流速1v 为2时渡河航线如下图：

0200400

60080010001200

-20

20

40

60

80

100

X: 29.99Y: 44.89

小船渡河航线

时间t

小船航线横、纵坐标值x 、y

X: 978.8

Y: -7.398e-07

X: 978.8Y: 50.03

从上图中可以看出，在流速1v 为2时, 当t=978.8秒时, 小船已到达对岸, 但是并没有到达B 点, 而是在B 点下游50.3米处。因为船头指向B 点, 即船头指向逆流方向, 且船在静水中的速度等于水流速度, 所以船永远到达不了B 点，即渡河时间为无穷大。

五、模型分析

用3/2阶龙格-库塔方法求得的小船航线的曲线图很逼真、很符合实际情况。随着1v 的变化，小船的航线也在变化，具体如下： 1). 当1v =0时，航线为直线；

2). 当0<1v <2时，航线呈一个类似于抛物线的曲线，1v 越大，“类抛物线”的顶点的横坐标越大，纵坐标也越大。即水速越大，航线的顶点顺流而下的距离就越大；

3). 当1v =2时，船会静止在B 点下游近似d /2的地方； 上述出现的情况都与我们的直观感觉相符合。