湖水污染问题的数学建模与求解

中国传媒大学2010 学年第一学期数学建模与数学实验课程

数学建模与数学实验

题目Pristine湖污染问题的建模与求解

学生姓名

学号

班级

学生所属学院

任课教师

教师所属学院

成绩

Pristine湖污染问题的建模与求解

摘要

本文讨论了湖水污染浓度变化趋势的预测问题。

通过分析水流输入输出湖泊的过程，建立了湖水污染浓度随时间变化的含参变量的微分方程模型，在河水污染浓度恒定和自然净化速率呈线性关系的情况下，求得其精确解，带入具体数据得到结论：在PCA声称的河水污染浓度下，湖的环境不会恶化；在工作人员实地测得的河水浓度下，湖的环境将会恶化。

同时建立了计算机模拟模型，带入具体数值，运用时间步长法来仿真模拟了在湖水污染浓度稳定以前湖水每天的变化情况，输出自PCA建厂以来每年的湖水污染浓度，得到与微分方程模型相同的结论。

在全停产和半停产时，通过前面的两个模型可以计算湖水污染浓度在自然净化影响下的恢复到净化指标所需的年限。并可得到结论：在半停产状态下，在选定的自然净化速率常数的约束下，只有当河水污染浓度降至原来的3.15%（自然净化速率呈线性关系），4.7%（自然净化速率呈指数关系），才有可能使河水在100年内恢复至0.001mol/l，然后给出整改建议。

一、问题重述

Pure河是流入Pristine湖的唯一河流。50年前PCA公司在此河旁建起一个生产设施并投入运行。PCA将为处理的湖水排入河中，导致Pristine湖被污染。PCA公司声称：已排放的废水的标准多年从未改变切不会对湖的环境有影响。

10L，流入（流出）的水流速度为149.1L/年。

现已知：Pristine湖的湖容量为15

PCA公司声称河水污染浓度仅为0.001mol/L，自工厂以来没有改变过。

讨论下列问题：

（1）建立数学模型用PCA提供的公开数据判断湖的环境是否会恶化；

（2）以目前湖水污染浓度0.03mol/L，和河水污染浓度0.05mol/L为新数据判断湖的环境是否会恶化；

二、模型的合理假设和符号系统

2.1 模型的合理假设

（1）降水量和增发量相等；

（2）湖中流入量和流出量相等且一直未变；

（3）污水量远小于河水注入量，且污水与河水混合均匀；

（4）湖水混合均匀，且流入污水的扩散速度无限大；

（5）湖内除Pure河外，无其他污染源；

2.2 符号系统

0ρ：河水污染浓度mol/L ；

ρ：湖水污染物浓度mol/L ；

V ：湖泊容量15

10L ；

c ：自然净化速率mol/（L 。年）

μ：流入（流出）的水流速度149.1L/年；

t ：从PCA 建厂至考察时刻的时间段。

三、问题的分析

3.1 问题分析：

对于问题中几个词语的理解：

1.是否会恶化——湖的环境恶化即湖水污染浓度大于0.001mol/L ，要判断其是否会恶化，则需计算在某一污染物积累速度（分析影响此速度的因素）下，湖水能达到的最大污染浓度和其变化趋势，以及湖水经几年超过0.001mol/L ，经过几年达到最大污染浓度。

2.自然净化——自然净化是独立的生态系统进行自我调节的方式之一，是在空气，阳光，水和细菌的参与下，进行包括物理沉降，化学反应和生物转化三大方面的活动，其最终作用是将污染物转化为无害物质，从而净化生态系统。当河水输入湖泊并均匀混合之后，影响湖水污染浓度的唯一因素便是自然净化速度。 湖水污染问题水流的动态流程图：

此问题中，我们考察对象是湖水污染浓度的变化趋势：

1.在整改之前，其增加的趋势，超过净化指标0.001ml/L （即湖水恶化）的可能性和时限；

2.在全停产（无污染物输入）和半停产的情况下，其降低的趋势，达到净化指标的时限。

3.在整改之后，其增加的趋势，未定与净化指标之下某一水平的时限、

在前假设条件的基础之上，湖水容量不变，出河水外无其他的污染源，故我们可将湖泊作为一个封闭的生态系统，其简化的湖水被污染的动态过程为：受污河输入湖泊，河水与湖水均匀混合，受污河水进行自我的净化，湖水数出湖泊。湖中污染物的量直接决定了湖水污染浓度，而污染物的量受到以下两方面因素的影响：1.河水的污染浓度及其流入速度（根据已知此速度不变），2.湖水的自然净化速度，前者使其增加，后者使其减少（负增加）。问题一、二、三的实质都是要分析污染浓度的变化趋势，其去表便在于前一因素的不同。问题一中，河水污染浓度不变，恒为0.001mol/L ；问题二中，河水污染浓度可能会变化，受PCA 效益的影响而按一定规律波动；问题三中，在全停产或半停产的情况下，和硕污染浓度为0或减为问题二中的一部分。后一因素（自然净化速度）在三个问题中的作用都是相同的。

根据微积分的知识可知，在适当短的时间段之内，通过建立微分方程，可以将连续的过程离散化，从而可以得到湖水污染浓度与时间之间的关系式。

利用时间步长发，缩小步长值（从年到月到天），并与微分方程所得的精确解做出比较。

四、模型建立与求解

问题一：

根据PCA 的公开申明和所提供数据，可认为：河水污染浓度恒为0.001mol/L 。从存在自然净化和不存在自然净化两个方面考虑： （1）．在不考虑自然净化的情况下：

由于假设湖中流入量和流出量相等，而在经过与湖水均匀混合后，流出湖水污染浓度明显减小，故流出污染物的量小于流出污染物的量，污染物将在湖中沉积，从而使湖水污染浓度增加，当其增加至于输入的河水污染浓度相等时，河水污染浓度达到最大，并稳定在这一数值，不在增加。

建立湖水污染浓度随时间变化的微分方程模型：设在极短时间dt ，湖水污染浓度

增加d ρ，在将湖水被污染这一连续动态过程简化为离散的瞬间静止状态（如问题分析中所述）之后，根据湖中剩余量=输入量—输出量，我们可以列微分方程如下：

0dt V

d V dt ρμρρρ

μ??+?=++?

化简可得：0()

d dt V ρμ

ρρ=?- (1)

带入1510V =L ，141.910μ=?L/年和00.001ρ=mol/L 的数据，我们可得到ρ和t 的关系式如下：

0.190.19

0.001(1)t t e e ρ-??=??-

通过此关系式我们可知，当t=∞时，湖水污染浓度将趋近与0.001mol/L，即湖的环境不会恶化。

建立计算机模拟模型：湖水污染浓度的变化时有湖中污染物随时间的积累而引起的，这个逐步积累的过程我们可以用计算机进行仿真模拟，其实质为完成一个循环累加的过程，并可改变时间步长，如一年一年的累积，一月一月的累积，一天一天的累积，从而使我们的模拟值逐步精确，可与微分方程求得精确解比较，分析误差。为提高模拟结果的精确性和运算的效率，我们采用了逐天累加，数出年污染浓度的方式。

模拟程序见附件一，从后面的分析中我们可知：在此河水污染浓度恒定，无自然净化的最简单的情况下建立的模型是以后问题的基础，后面的问题只是改变条件或数据，其实质是不变的，故我们在此程序中加入了多个选择语句，在不同德条件或数据下执行不同的命令，从而用一个程序解决全部的模拟问题。

模拟所得的数据如下表1：

年份 1 3 5 8 10

ρ4

10- 1.73 4.35 6.14 7.82 8.51 年份15 20 30 40 50

ρ4

10-9.43 9.78 9.97 9.995 9.999

根据模拟数据所作的湖水污染浓度变化趋势图如下图1：

模拟所得的数据显示：湖水污染 浓度将稳定于0.001mol/L ，这与分析微分方程所得的结果是相符的。

（2）在考虑自然净化的情况下： A.首先，我们应该了解什么事自然净化。

根据资料显示：湖水中的污染物可分为有机污染物和无机污染物两大类，在多种环境因素（阳光，空气，水，水中生物，水中化学物质，重力等）的作用下，通过物理沉降，化学反应和生物转化一系列复杂的活动，它们的量会发生变化。 有机物在水环境中的迁移转化过程如图2：

可以看出整个过程是相当复杂的，不仅过程多，而且在相同的过程中，不同的物质有着不同的结果。此问题中没有明确给出输入湖中的污染物种类，也没有对湖泊环境做任何描述，无疑给问题的解决增加了极大的难度。为是问题简化，我们从一般的情况出发，假设：湖中所进行的反应均为一级反应，有机物的存在不会对环境参数造成改变，环境固定（自然净化速率恒定），考虑湖泊生态系统中起主要作用的几种过程作简要分析：

1.物理沉降。不同物质有着不同的沉降速度常数λ，沉降速率为c λρ=?；

2.挥发。不同物质有着不同的挥发速率常数v

K ，在有机物在水体上的大气中的

分压为0的条件下，挥发速率为

v K c Z ρ

?=

(Z 为水体深度)

3.水解反应。不同物质有着不同的水解速率常数

p

K ，在一级反应的条件下，水

解速率为

p c K ρ

=?；

4.生物降解反应。不同物质有着不同的降解速率常数b

K ，在一级反应的条件下，

生物降解速率

b c K ρ

=?；

B.在最简模拟程序的基础上，从每天的积累量中减去每天的自然净化量K ?时间

步长，重复循环可得到逐年湖水污染浓度的值，设1

0.0003K -=天。K 越大，湖

水污染浓度将越快稳定于一个更小的浓度值。 模拟所得的数据如下表2：

年份

1 3 5 8 10 ρ 410-

1.587 3.437 4.315 4.878 4.984 年份

15 20 30 40 50 ρ 410-

5.087

5.103

5.106

5.106

5.106

根据模拟数据所作的湖水污染浓度变化趋势图如下图3：

C.进一步考虑自然净化速率的影响。更贴近于实际的情况是，自然净化速率c 与湖水污染浓度ρ成指数关系，c 随着ρ的增加而增加，但增加的速率会逐渐减小，

用关系式表达即可为：

(1)B c A e ρ

?=-，A,B 为与环境和污染物种类有关的常数。将此关系式带入微分方程（2），得到一个新的微分方程，此方程无解。但是我们

可以通过计算机模拟求得数值解。

A 的大小影响着ρ最终稳定浓度，

B 的大小影响着ρ达到最终稳定浓度的快慢。在前面达到线性关系的基础上，在ρ=0/。001mol/L 处的c 值大小

确定A 至数量级，并稳定年限定在10—20年间，从而得到一组估计值：

4A 810,B 400-=?=

模拟所得数据如下表3：

年份

1 3 5 8 10 ρ 410-

1.493

2.971

3.539 3.816 3.871 年份

15 20 30 40 50 ρ 410-

3.903

3.906

3.906

3.906

3.906

根据模拟所得数据所作的湖水污染浓度变化趋势图如下图4：

观察数据可知：无论自然净化速率c 与污染浓度ρ成线性关系还是指数关系，当河水污染浓度

0ρ=0.001mol/L 时，ρ都将稳定于一个小于0.001mol/L 的值，

也即：

湖的环境不会恶化。 问题二：

根据实际情况，在此我们只考虑存在自然净化的情况。

（1）.河水污染浓度恒定为0.05mol/L: 与问题一（2）相同，只是

0ρ=0.05mol/L ，沿用问题一（2）的微分方程模型和差

分模拟模型，我们可以得到以下结果：

A.微分方程模型（自然净化速率c 与湖水污染浓度ρ成线性关系）：

当时间t =∞时，

0KV ρμ

ρμ=

+，即湖水污染浓度稳定于一个与K 有关的值。当

0ρ=0.05mol/L 时，如要ρ稳定于0.001mol/L 则K 值为1-9.31年。根据我们设定的1

K 0.0003-=天可计算，湖水最终污染浓度

0ρ=0.0317mol/L ，超过净化指标

0.001mol/L ，故在此条件下，湖的环境将会恶化。 B.差分模拟模型：

在c 与ρ成线性关系时，得到模拟数据如下表4：

年份

1 3 5 8 10 ρ 210-

0.82 1.88 2.46 2.88 3.01 年份

15 20 30 40 50 ρ 210-

3.13

3.16

3.17

3.17

3.17

在c 与ρ成指数关系时，得到模拟数据如下表5：

年份

1 3 5 8 10 ρ 210-

0.81 2.01 2.81 3.59 4.02 年份

15 20 30 40 50 ρ 210-

4.32

4.48

4.56

4.58

4.58

观察以上数据可知：由于河水污染没得过大，河水污染没得从第一年起就超过了净化指标，并逐年增加，是湖的环境恶化。 （2）.河水污染浓度变化：

根据实际情况，一个工厂的生产量并非恒定不变，每年每月甚至每天也有所不同，从而起其排放污染物的量也将随生产量的变化而变化的。假设PCA 自建立以来的年排污量服从Logistic 模型，并考虑到污染物的量受到多方面因素的影响，在

每一时刻的梁上加上一个服从正太分布，范围在此量的10%内的较小量ξ，则可建立河水污染浓度的模型如下：

01r M t

M

Ce ρξ??=

++（M,C,r 为与工厂效益有关的常数）

设此时工厂处于文本上升的发展阶段，其变化曲线如下图5：

通过计算机模拟产生随机数（如图中星点所示），带入模拟程序。 在c 于ρ成线性关系时，得到的模拟数据如下表6： 年份

1 5 10 20 2

2 35 40 ρ

410-

0.185

0.801

1.792

7.272

9.581

52.04

97.13

年份 50 82 90 100 112 140 150 ρ

410-

260.7

625.0

633.9

617.2

626.3

628.5

632.5

根据模拟所得数据所作的湖水污染浓度变化趋势如下图5：

在c 于ρ成指数关系时，得到的模拟数据如下表7： 年份

1 5

10 20 24 35 45 ρ

410-

0.169

0.586

1.250

5.194

9.046

44.26

95.75

年份 50 80 100 120 130 140 150 ρ

410-

321.7

925.5

959.6

950.6

949.6

952.5

964.8

根据模拟所得数据所作的湖水污染浓度变化趋势如下图5：

在第50年时，无论自然净化速率c 于湖水污染浓度ρ成线性关系还是指数关系，可以看到湖水污染浓度都接近工作人员实地测得值，这从检验的角度说明我们的对模型参数的估计是可取的。

观察模拟所得数据可知，湖的环境将会恶化。

五、模型的评价及评改进

5.1模型的评价

（1）.对于首先建立的含参数的微分方程模型，在其有精确解的时候，能给出湖水污染浓度随时间变化的规律，得到很好的预测效果，但微分方程并不是任何条件下都有解，在本问题中，只有当河水污染浓度恒定且自然净化速率呈线性关系时，微分方程可解，故适用的范围受到了限制。 （2）.差分模拟法，即用数值方法求解微分方程，它解决了微分方程无精确解的 情况，适用范围很大。使用一个兼顾运算效率和模拟结果精度的时间步长是该方 法的关键。在处理本问题时，采用的时间步长为天。这样处理的精度大约为0.0001, 与精确解相当接近。

5.2模型的改进

（1）.自然净化速率的改进：资料显示，污染物的很多生物降解过程是通过好氧 反应实现的，当污染浓度较大时，需氧量增加，而湖中氧气总量不变，自然净化 速率会有所减小，所以，自然净化速率的函数可作为分段函数考虑。前一段还是

用前面的指数模型，后面一段可以考虑为

(,,0)D E

C e C

D

E ρ-?+?>。模型这样改进后，更加符合实际情况。

（2）.河水浓度变化规律模型的改进：根据价值规律，一个工厂的发展总是有起 有落，其生产量绕某一中心值上下波动，故生产河水的浓度的变化应该是具有一 定周期性，同时又有所增加。可以考虑在原来的Logistic 模型后加上一个傅立叶 函数（其周期就是工厂生产的周期）和一个正态随机变量。

六、参考文献

[1]《数学建模与数学实验》 高等教育出版社，赵静、但琦主编，2008年1月第三版

[2]《环境化学》，高等教育出版社出版，戴树桂主编，1997年3月；

七、附录

附录一 :

问题一：

（1）河水浓度恒为0.001mol/L ，不考虑自然净化，

a：采用差分方程模型

命令：wenti12（‘d’，0.50,0.001，‘n’）

b：采用微分方程模型

命令：t=1:1:51；

P=0.001.*exp(-0.19*t)

P

plot(t,p)

grid on

（2）河水浓度恒为0.001mol/L，考虑自然净化，

a：采用差分自然净化

线性自然净化

命令：wenti12('d',0,50,0.001,'y','xx')

指数自然净化

命令：wenti12('d',0,50,0.001,'y','zs')

b：采用微分方程模型，线性自然净化

命令：t=1:1:51;

p=0.001.*exp(-0.19.*t).*(-1+1.*exp(0.19.*t));

p

hold on

plot(t,p)

grid on

问题二：

（1）河水浓度恒为0.05mol/L,不考虑自然净化

a：采用差分方程模型

命令：wenti12('d',0.50,0.05,'n')

b:采用微分方程模型

命令：t=1:1:51;

p=0.05.*exp(-0.19.*t).*(-1+1.*exp(0.19.*t));

p

hold on

plot(t,p)

grid on

（2）河水浓度恒为0.05mol/L，考虑自然净化，

a:采用差分方程模型

线性自然净化

命令：wenti12('d',0.50,0.05,'y','xx')

指数自然净化

命令：wneti12('d',0.50,0.05,'y','zs')

b：采用微分方程模型，线性自然净化

命令

t=1:1:51;

p=0.00166945.*exp(-0.2995.*t).*(-19+19.*exp(0.2995.*t));

p

hold on

plot(t,p)

grid on

（3）河水浓度变化，

不考虑自然净化

命令：wenti12('d',0.50,-1,'n')

考虑自然净化

线性自然净化

命令：wenti12('d',0,50,-1,'y','xx')

指数自然净化

命令：wenti12('d',0,50,-1,'y','zs')

附录二：

主程序：

function wenti12(l,chu,n,s,w,b)

%按照时间步长法，以每年或每月或每日为时间段，求出浓度的变化，最后输出每年的浓度，做一个大致的观察,

%格式：wenti12（P1,P2,P3,P4,P5,P6)

%P1:y,m,d,h,min,s 表示设定步长值

%P2：湖水污染初值

%P3：年份

%P4：河流污染初值

%P5: 是否考虑自然降解,n-不考虑,y-考虑

%P6: xx 采用线性净化模型，zs 采用指数净化模型

format compact

format long

if nargin 2 %50年,河水污染变化,不考虑自然降解

n=50;

s=-1;

w='n';

b='xx';

elseif nargin 3 %河水污染变化,不考虑自然降解

s=-1;

w='n';

b='xx';

end

if l 'y'

q=1

elseif l 'm'

q=12

elseif l 'd'

q=365

elseif l 'h'

q=365*24

elseif l 'min'

q=365*24*60

elseif l 's'

q=365*24*60*60

else

error('error,wrong parameter...try again!')

end

A=zeros(1,n+1);

B=zeros(1,n+1);

A(1)=chu;

B(1)=chu;

var=0; %自然降解

var1=0; %河水污染

if w 'n'

var=0;

elseif w 'y'

var=fen(A(1),l,b);

else

error('error,wrong parameter..try again!')

end

for i=1:1:n

if s -1

var1=river(i);

else

var1=s;

end

for j=1:1:q

if w 'n'

var=0;

elseif w 'y'

var=fen(A(i),l,b);

end

A(i)=(A(i)*10^15+var1*1.9*10^14/q)/(10^15+1.9*10^14/q)-var;

B(i)=0.19*var1/q+(1-0.19/q)*B(i)-var;

end

if w 'n'

var=0;

elseif w 'y'

var=fen(A(i),l,b);

end

A(i+1)=(A(i)*10^15+var1*1.9*10^14/q)/(10^15+1.9*10^14/q)-var;

B(i+1)=0.19*var1/q+(1-0.19/q)*B(i)-var;

end

[A;B]

hold on

plot(1:1:n+1,A,'o')

grid on

function p=wenti3(l,bi,year,b,a)

%格式: wenti3(P1,P2,P3,P4，P5)

%P1:y,m,d,h,min,s 表示设定步长值

%P2:百分比

%P3:可认为的最大净化年限

%P4:xx-线性自然净化,zs-指数自然净化

%P5：湖水浓度初始值

format compact

format long

if l 'y'

q=1

elseif l 'm'

q=12

elseif l 'd'

q=365;

elseif l 'h'

q=365*24

elseif l 'min'

q=365*24*60

elseif l 's'

q=365*24*60*60

else

error('error,try again!')

end

A=zeros(100);

B=zeros(100);

temp=0.05*bi;

A(1)=a;

for i=1:1:year

for j=1:1:q

A(i)=(A(i)*10^15+temp*1.9*10^14/q)/(10^15+1.9*10^14/q)-fen(A(i),l,b);

t=A(i);

if A(i)<0.001 break;

end

end

A(i+1)=(A(i)*10^15+temp*1.9*10^14/q)/(10^15+1.9*10^14/q)-fen(A(i),l,b );

if A(i)<0.001 break;

end

end

if i year

p=0;

else

p=1;

end

%程序调试所用

%sprintf('停产')

sprintf('需要%d 年%d天~~',i,j)

%sprintf('第%d年%d天为:%.12f',i,j-1,t)

%sprintf('第%d年%d天为:%.12f',i,j,A(i))

function low=wenti3_2(b,year,l,h)

%采用二分法求出百分比

format compact

format long

low=l;

high=h;

cen=(low+high)/2

cha=high-low;

while cha>10^-6 %误差判断条件

m=wenti3('d',cen,year,b,0.05); %调用函数

if m 0

high=cen;

cen=(low+high)/2

else

low=cen;

cen=(low+high)/2

end

cha=high-low;

end

high;

low;

sprintf('工厂在生产为原来的%f时，湖水才能在%d年内净化。',low,year) function wenti3_3(l,year,bi,b)

%格式: wenti3_3(P1,P2,P3,P4)

%P1:y,m,d,h,min,s 表示设定步长值

%P2:year

%P3: 百分比

%P4:xx-线性自然净化,zs-指数自然净化

format compact

format long

if l 'y'

q=1

elseif l 'm'

q=12

elseif l 'd'

q=365;

elseif l 'h'

q=365*24

elseif l 'min'

q=365*24*60

elseif l 's'

q=365*24*60*60

else

error('error,try again!')

end

A=zeros(1,year);

A(1)=0.03;

N=0;

for i=1:1:year

temp=0.05*bi;

for j=1:1:fix(q)

A(i)=(A(i)*10^15+temp*1.9*10^14/q)/(10^15+1.9*10^14/q)-fen(A(i),l,b);

end

A(i+1)=(A(i)*10^15+temp*1.9*10^14/q)/(10^15+1.9*10^14/q)-fen(A(i),l,b );

end

N=A(i)

被调用的程序：

function y=fen(N,l,b)

%自然净化

%格式:fen(P1,P2)

%P1:前一时间点的污染值

%P2:y,m,d,按年或月或日

%线性关系

if b 'xx'

first=0.1095;

%first=0.931;

if l 'y'

y=first*N;

elseif l 'm'

y=first*N/12;

elseif l 'd'

y=first*N/365;

else

error('error,wrong parameter...try again!') end

elseif b 'zs'

%指数关系

second=8*10^-4;

if l 'y'

y=second*(1-exp(-N*400));

elseif l 'm'

y=second*(1-exp(-N*400))/12;

elseif l 'd'

y=second*(1-exp(-N*400))/365;

else

error('error,wrong parameter...try again!') end

else

error('error,wrong parameter...try again!') end

%河水污染浓度函数的参数方程，

function F=fun(x)

M=x(1);

C=x(2);

k=x(3);

F=[eps-M/(1+C);

0.05-M/(1+C*exp(-k*M*50));

0.1- M/(1+C*exp(-k*M*100))];

function y=river(n)

%河水浓度变化

%本函数是求第n年的河水污染浓度,返回值是第n年的河水污染浓度

%分为三个阶段；

%一：%河水污染浓度是一个关于时间n的函数，先求出此函数的三个未知参数M，C，k

%二：%在所求出的函数上加上一个小的波动，更能满足实际要求

%三：%计算第n年的河水污染浓度

%一：求参数

%由于要提供搜索点，所以先在mathmatica中用solve命令求出大致的M ，C，k 值，作为搜索点

%mathematica中的程序为：

% Solve[{0.0001 M/(1 + C),

% 0.05 M/(1 + C*E^(-k*M*50)),

% 0.1 M/(1 + C*E^(-k*M*100))}, {M, C, k}]

%x0=[0.1001;1000;1.37997]; %搜索点

%options=optimset('Display','iter');

%[x,fval]=fsolve('fun',x0,options) %解方程组

x=1000*[0.00010010030107 1.00000021746742 0.00137976638106];

i=1:1:100;

g=x(1)./(1+x(2).*exp(-x(1).*x(3).*i));

%plot(i,g)

%二:加上一个波动,10%的波动范围

h=zeros(1,100);

sigma=0.03;

for j=1:100

para=sigma*randn;

h(j)=g(j)*(1+para);

end

j=1:100;

%画出波动图,调试时用

%plot(i,g,j,h,'xr')

%计算第n年的浓度

y=x(1)./(1+x(2).*exp(-x(1).*x(3).*n));

y=y*(1+sigma*randn);

数学建模实验答案-概率模型

数学建模实验答案-概率模型

实验10 概率模型（2学时） （第9章 概率模型） 1.（验证）报童的诀窍p302~304, 323（习题2） 关于每天报纸购进量的优化模型： 已知b 为每份报纸的购进价，a 为零售价，c 为退回价（a > b > c ），每天报纸的需求量为r 份的概率是f (r )(r =0,1,2,…)。 求每天购进量n 份，使日平均收入，即 1 ()[()()()]()()()n r r n G n a b r b c n r f r a b nf r ∞ ==+=----+ -∑∑ 达到最大。 视r 为连续变量，f (r )转化为概率密度函数p (r )，则所求n *满足 * ()n a b p r dr a c -= -? 已知b =, a =1, c =，r 服从均值μ=500（份），均方差σ=50（份）的正态分布。报童每天应购进多少份报纸才能使平均收入最高，这个最高收入是多少 [提示：normpdf, normcdf] 要求：

(1) 在同一图形窗口内绘制10 ()()n y n p r dr =?和2()a b y n a c -= -的图形，观察其交点。 [提示] 22 ()2()r p r μσ-- = ，0 ()()()n n p r dr p r dr p r dr -∞ -∞ =-?? ? ☆(1) 运行程序并给出结果： (2) 求方程0()n a b p r dr a c -= -?的根n *（四舍五入取整），并求G (n *)。

mu=500;sigma=50; a=1; b=; c=; r=n+1; while (a-b)*n*normpdf(r,mu,sigma)>1e-6 r=r+1; end r=n+1:r; G=sum((a-b)*n*normpdf(r,mu,sigma)); r=0:n; G=G+sum(((a-b)*r-(b-c)*(n-r)).*normpdf(r,mu,sigma)) ☆(2) 运行程序并给出结果： 2.（编程）轧钢中的浪费p307~310 设要轧制长l=的成品钢材，由粗轧设备等因素决定的粗轧冷却后钢材长度的均方差σ=，问这时钢材长度的均值m应调整到多少使浪费最少。 平均每得到一根成品材所需钢材的长度为 () () m J m P m = 其中， 2 2 () 2 ()(), () 2 x m l P m p x dx p xσ πσ - - ∞ == ? 求m使J(m)达到最小。 等价于求方程 () () z z z λ ? Φ =- 的根z*。 其中：

空气污染数学建模.docx

A.污染气体的传播扩散 摘要 钢铁生产排放的污染气体是造成雾霾的重要原因之一，研究污染气体的扩散特征，正确模拟污染气体的扩散过程，能够为钢铁生产集团提出更好的治理管理措施，具有实际意义。 针对问题一：污染气体的排放速度为300m/s,在不考虑风向风速及高度影响的情况下，此问题即为二维平面的连续点源扩散问题，由此在二维xoy 平面上建立连续点源扩散方程模型，其中为气体扩散系数，本文中取为常数10，f(x,y,t ) 为污染气体的排放速度，在本文中恒为300m/s ；对上述偏微分方程模型，本文采用ADI法（Alternating direction implicit，交替方向隐式法）求解出迭代格式，利用MATLAB编程，求出模型一的数值解，并得到任意时刻污染气体的浓度分布情况。通过SPSS软件，对附件一所给的原始实际数据与模型一求解得到的模拟值进行显着性检验，检验结果显示该模型与实际情况吻合。 针对问题二：考虑风向风速对污染气体扩散过程的影响时，在基于对问题一求解的基础上，在模型一的扩散方程模型中加入风向风速的平流项，由此得到有风情况下的模型，其中分别为风速在x, y方向的分量；对此模型同样采用ADI法求出迭代格式，利用MATLAB编程，求出模型二的数值解，并得到任意时刻污染气体的浓度分布情况。通过SPSS软件，对附件二所给的原始实际数据与模型二求解得到的模拟值进行显着性检验，检验结果显示该模型与实际情况吻合。 针对问题三：考虑有风时增加高度的影响，此问题即为三维空间的污染气体扩散问题，考虑到三维模型的编程复杂度，而且污染气体的扩散在xoy平面上各向同性，可以将污染气体在y方向的扩散等价为在x方向上的扩散，此时便只需要建立xoz平面上的扩散模型。在基于对问题二求解的基础上，在模型二的扩散方程中增加高度项，由此得到模型三为，其中为z 方向的扩散系数；对该扩散方程同样采用ADI法求出迭代格式，利用MATLAB编程，求出模型二的数值解，并得到任意时刻污染气体的浓度分布情况。关键词：污染气体扩散方程ADI法数值解 一、问题重述 目前，治理雾霾是人们最为关心的热点问题之一。中国社科院发布的《气候变化绿皮书》中提及，雾霾形成的原因里，重工业、车辆尾气、土方施工都榜上有名，其中钢铁生产也是造成雾霾的重要原因之一。 某钢铁生产集团烟囱污染气体的排放对周边地区大气污染的影响非常大，为了提出更好的治理管理措施，需要对其污染气体扩散的特征进行分析。现在，我们需要在三种情况下考虑污染气体的扩散过程： 1.在不考虑风向和高度影响的情况下，建立模型，模拟某钢铁生产集团的烟囱排放污染 气体的扩散过程，假设烟囱的排放速度为300m/s。 2.考虑风向为东北风，平均风速0.6m/s的情况下，模拟污染气体的传播扩散过程。 3.在考虑风向的基础上增加高度的影响，建立模型，模拟污染气体的传播扩散过程。 4.基于上述模型结论，给该钢铁生产集团提供一个污染气体治理建议报告。

《数学建模与数学实验》本科教学日历

《数学建模与数学实验》本科教学日历 数学建模部分 开设课程课程名称数学建模课程编号0701107 施教单位理学院 课内学时 总课时36 课程性质公共基础讲授课时28 修读要求选修实践课时8 选用教材教材名称数学建模教程出版社名称高等教育出版社 出版时间 及版次 2011年出版，第一版印刷时间2011年 其他情况 教学安排 班次授课对象及人数任教教员（指导教员）姓名及职称数学建模A 各专业本科学员 吴孟达教授 段晓君教授 毛紫阳讲师 王丹讲师 数学建模B 各专业本科学员 吴孟达教授 段晓君教授 毛紫阳讲师 王丹讲师 课次节 次 授课内容 教学 方法 采用现代化教学手段(课时) 多媒体电教双语网络实验 1 1 （1）什么是数学建模？数学建模的一般概念 （2）几个数学建模问题 讲授 1 2 （1）数学建模的一般步骤 （2）敏感问题调查案例 讲授 1 2 3 （1）行走步长问题 （2）雨中行走淋雨量最小问题 （3）道路是越多越通畅吗? 讲授 1 4 （1）有奖销售的抽奖策略问题 （2）“非诚勿扰”女生最佳选择问题 （3）网络文章流行度预测和招聘匹配 讲授 1 3 5 （1）线性规划模型基本概念 （2）整数规划模型 （3）0-1规划模型 讲授 1 6 （1）非线性规划 （2）多目标规划 讲授 1 4 7 （1）最短路算法 （2）最小生成树算法 讲授 1 8 （1）最大流算法 （2）PageRank算法 讲授 1 5 9 规划模型上机实践实践 1

课次节 次 授课内容 教学 方法 采用现代化教学手段(课时) 多媒体电教双语网络实验10 图论模型上机实践实践 1 6 11 （1）博弈模型基本概念 （2）Nash平衡和Pareto最优 （3）博弈论案例 讲授 1 12 （1）贝叶斯纳什均衡 （2）拍卖模型 讲授 1 7 13 社会选择理论中的选举问题数学模型-阿罗不可能定理讲授 1 14 越野长袍团体赛排名规则公平性问题讲授 1 8 15 军事作战模型-Lanchester作战模型讲授 1 16 自动化车床管理模型讲授 1 9 17 （1）“边际效应”基本概念 （2）实物交换模型，最佳消费模型、报童售报问题 讲授 1 18 （1）价格弹性模型 （2）合作效益的Shapley值分配模型 讲授 1 10 19 （1）聚类分析基本概念 （2）常用聚类算法 讲授 1 20 （1）方差分析基本概念 （2）单因素方差分析 （3）双因素方差分析 讲授 1 11 21 （1）主成分分析基本概念 （2）因子分析 讲授 1 22 （1）一元回归分析 （2）多元回归分析 （3）多元回归模型的检验与优化 讲授 1 12 23 聚类分析和方差分析上机实践实践 1 24 主成分分析和多元回归分析上机实践实践 1 13 25 （1）遗传算法基本思想 （2）算法步骤 讲授 1 26 遗传算法计算实例讲授 1 14 27 （1）模拟退火算法基本思想 （2）算法步骤 讲授 1 28 模拟退火算法计算实例讲授 1 15 29 （1）蚁群算法基本思想 （2）算法步骤 讲授 1 30 （1）数学建模中的计算机仿真 （2）不可召回的秘书招聘问题 （3）车灯光源优化设计 （4）生命游戏 讲授 1 16 31 遗传算法上机实践实践 1 32 模拟退火算法上机实践实践 1

数学建模与数学实验习题

数学建模与数学实验课程总结与练习内容总结 第一章 1.简述数学建模的一般步骤。 2.简述数学建模的分类方法。 3.简述数学模型与建模过程的特点。 第二章 4.抢渡长江模型的前3问。 5.补充的输油管道优化设计。 6.非线性方程（组）求近似根方法。 第三章 7.层次结构模型的构造。 8.成对比较矩阵的一致性分析。 第五章 9.曲线拟合法与最小二乘法。 10 分段插值法。 第六章 11 指数模型及LOGISTIC模型的求解与性质。 12.VOLTERRA模型在相平面上求解及周期平均值。 13 差分方程（组）的平衡点及稳定性。 14 一阶差分方程求解。 15 养老保险模型。

16 金融公司支付基金的流动。 17 LESLLIE 模型。 18 泛函极值的欧拉方法。 19 最短路问题的邻接矩阵。 20 最优化问题的一般数学描述。 21 马尔科夫过程的平衡点。 22 零件的预防性更换。 练习集锦 1. 在层次分析法建模中，我们介绍了成对比较矩阵概念，已知矩阵P 是成对比较矩阵 31/52a b P c d e f ?? ??=?????? ，（1）确定矩阵P 的未知元素。 （2）求 P 模最大特征值。 （3）分析矩阵P 的一致性是否可以接受（随机一致性指标ＲＩ取0.58）。 2. 在层次分析法建模中，我们介绍了成对比较矩阵概念，已知矩阵P 是三阶成对比较矩阵 322P ? ???=?????? ，（1）将矩阵P 元素补全。 （2）求P 模最 大特征值。 （3）分析矩阵P 的一致性是否可以接受。 3.考虑下表数据

(1)用曲改直的思想确定经验公式形式。 （2）用最小二乘法确定经验公式系数。 4.. 考虑微分方程 (0.2)0.0001(0.4)0.00001dx x xy dt dy y xy dt εε?=--????=-++?? （1）在像平面上解此微分方程组。（2）计算0ε=时的周期平均值。（3）计算0.1ε=时，y 的周期平均值占总量的周期平均值的比例增加了多少？ 5考虑种群增长模型 '()(1/1000),(0)200x t kx x x =-= （1）求种群量增长最快的时刻。（2）根据下表数据估计参数k 值。 6. 布均匀，若环保部门及时发现并从某时刻起切断污染源，并更新湖水（此处更新指用新鲜水替换污染水），设湖水更新速率是 3 (m r s 单位：）。 （1） 试建立湖中污染物浓度随时间下降的数学模型？ 求出污染物浓度降为控制前的5%所需要的时间。 7. 假如保险公司请你帮他们设计一个险种:35岁起保,每月交费400元,60岁开始领取养老金,每月养老金标准为３６００元，请估算该保险费月利率为多少（保留到小数点后５位）？ 8. 某校共有学生40000人，平时均在学生食堂就餐。该校共有,,A B C 3 个学生食堂。经过近一年的统计观测发现：A 食堂分别有10%，25%的学生经常去B ，C 食堂就餐，B 食堂经常分别有15%，25%的同学去

2011年数学建模-重金属污染模型(特等奖作品)

兰州理工大学数学建模面朝大海 2011高教社杯全国大学生数学建模竞赛 承诺书 我们仔细阅读了中国大学生数学建模竞赛的竞赛规则. 我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮 件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问 题。 我们知道，抄袭别人的成果是违反竞赛规则的, 如果引用别人的成果或其他 公开的资料（包括网上查到的资料），必须按照规定的参考文献的表述方式在正 文引用处和参考文献中明确列出。 我们郑重承诺，严格遵守竞赛规则，以保证竞赛的公正、公平性。如有违反 竞赛规则的行为，我们将受到严肃处理。 我们参赛选择的题号是（从A/B/C/D中选择一项填写）： 我们的参赛报名号为（如果赛区设置报名号的话）： 所属学校（请填写完整的全名）：兰州理工大学 参赛队员 (打印并签名) ：1. 2. 3. 指导教师或指导教师组负责人 (打印并签名)： 日期： 2011 年 9 月 11 日赛区评阅编号（由赛区组委会评阅前进行编号）：

2011高教社杯全国大学生数学建模竞赛 编号专用页 赛区评阅编号（由赛区组委会评阅前进行编号）： 全国统一编号（由赛区组委会送交全国前编号）：全国评阅编号（由全国组委会评阅前进行编号）：

城市表层土壤重金属污染分析 摘要 随着城市经济的快速发展和城市人口的不断增加，人类活动对城市环境质量的影响日显突出。本文根据大气中重金属污染物的传播特征，将高斯扩散模型应用于土壤重金属污染问题，建立数学模型并进行求解。 问题一中，利用题目所提供的重金属元素样本值，借助matlab软件，形象地绘制出了城市地形图以及八种重金属在该城市五个功能区的二维等高线分布图。对所收集的土壤样品分别采用单因子污染指数法与内梅罗（N.L.Nemerow）综合污染指数法对该城市分区进行污染程度评价。 问题二中，通过对所给数据的分析，利用主成分分析法，推断出了重金属污染的主要原因。 问题三中，通过数据分析，发现了金属污染物的传播特征，由此建立重金属高斯扩散模型。利用数形结合的方法最终确定了污染源的位置。 问题四中，通过对建模过程中考虑因素的全面性以及所忽略因素的得当性分析，评估了所建模型的优缺点，确定了更好地研究城市地质环境的演变模式所需要搜集的信息，并给出了更合理建立数学模型的方法。 关键词：重金属污染主成分分析法单因子指数法内梅罗综合污染指数法高斯模型污染源

数学建模实验答案初等模型

实验02 初等模型（4学时） （第2章初等模型） 1.（编程）光盘的数据容量p23~27 表1 3种光盘的基本数据 CAV光盘：恒定角速度的光盘。 CLV光盘：恒定线速度的光盘。 R2=58 mm, R1=22.5 mm，d, ρ见表1。

CLV光盘的信息总长度(mm) L CLV 22 21 () R R d π- ≈ CLV光盘的信息容量(MB) C CLV = ρL CLV / (10^6) CLV光盘的影像时间(min) T CLV = C CLV / (0.62×60) CAV光盘的信息总长度(mm) L CAV 2 2 2 R d π≈ CAV光盘的信息容量(MB) C CAV = ρL CAV / (10^6) CAV光盘的影像时间(min ) T CAV = C CAV / (0.62×60) 1.1（验证、编程）模型求解 要求： ①（验证）分别计算出LCLV, CCLV和TCLV三个3行1列的列向量，仍后输出结果，并与P26的表2（教材）比较。 程序如下：

②（编程）对于LCAV, CCAV和TCAV，编写类似①的程序，并运行，结果与P26的表3（教材）比较。 ★要求①的程序的运行结果： ★要求②的程序及其运行结果：

1.2（编程）结果分析 信道长度LCLV 的精确计算：21 2R CLV R L d π=? 模型给出的是近似值：2221() CLV R R L L d π-= ≈ 相对误差为：CLV L L L δ-= 要求：

①取R2=58 mm, R1=22.5 mm，d, ρ见表1（题1）。 分别计算出LCLV, L和delta三个3行1列的列向量，仍后将它组合起来输出一个3行3列的结果。 ②结果与P26的表2和P27（教材）的结果比较。 [提示] 定积分计算用quad、quadl或trapz函数，注意要分别取d的元素来计算。要用数组d参与计算，可用quadv（用help查看其用法）。 ★编写的程序和运行结果： 程序：

数学建模与数学实验课后习题答案

P59 4.学校共1002名学生，237人住在A 宿舍，333人住在B 宿舍，432人住在C 宿舍。学生要组织一个10人的委员会，使用Q 值法分配各宿舍的委员数。 解:设P 表示人数，N 表示要分配的总席位数。i 表示各个宿舍（分别取A,B,C ），i p 表示i 宿舍现有住宿人数，i n 表示i 宿舍分配到的委员席位。 首先，我们先按比例分配委员席位。 A 宿舍为:A n = 365.21002 10237=? B 宿舍为:B n =323.31002 10333=? C 宿舍为:C n =311.4100210432=? 现已分完9人，剩1人用Q 值法分配。 5.93613 22372 =?=A Q 7.92404 33332 =?=B Q 2.93315 44322 =?=C Q 经比较可得，最后一席位应分给A 宿舍。 所以，总的席位分配应为：A 宿舍3个席位，B 宿舍3个席位，C 宿舍4个席位。

商人们怎样安全过河

由上题可求：4个商人，4个随从安全过河的方案。 解：用最多乘两人的船，无法安全过河。所以需要改乘最多三人乘坐的船。 如图所示，图中实线表示为从开始的岸边到河对岸，虚线表示从河对岸回来。商人只需要按照图中的步骤走，即可安全渡河。总共需要9步。

P60 液体在水平等直径的管内流动，设两点的压强差ΔP 与下列变量有关：管径d,ρ,v,l,μ,管壁粗糙度Δ，试求ΔP 的表达式 解：物理量之间的关系写为为()?=?,,,,,μρ?l v d p 。 各个物理量的量纲分别为 []32-=?MT L p ，[]L d =，[]M L 3-=ρ，[]1-=LT v ，[]L l =，[]11--=MT L μ，Δ是一个无量纲量。 ???? ??????-----=?0310100011110010021113173A 其中0=Ay 解得 ()T y 00012111---=， ()T y 00101102--=， ()T y 01003103--=， ()T y 10000004= 所以 l v d 2111---=ρπ，μρπ112--=v ，p v ?=--313ρπ，?=4π 因为()0,,,,,,=??p l v d f μρ与()0,,,4321=ππππF 是等价的，所以ΔP 的表达式为： ()213,ππψρv p =?

数学建模_湖水污染问题(1)

湖水污染问题 一.问题提出 下图是一个容量为2000m3的一个小湖的示意图，通过小河A水以 /s的速度流入，以相同的流量湖水通过B流出。在上午8:00，因交通事故，一辆运输车上一个盛有毒性化学物质的容器倾翻，在图中X点处注入湖中。在采取紧急措施后，于上午9：00事故得到控制，但数量不详的化学物质Z已泻入湖中，初步估计Z的数量在5m3至20m3之间。 （1）请建立一个数学模型，通过它来估计湖水污染程度随时间的变化； （2）估计湖水何时到达污染高峰； （3）何时污染程度可降至安全水平(<=%)。 二.模型假设 1、湖水流量为常量，湖水体积为常量； 2、流入流出湖水水污染浓度为常量 三.问题分析 分析：湖水在时间t时污染程度，可用污染度F（t）表示，即每立方米受污染的水中含有Fm3的化学污染物质和（1-F）m3的清洁水。用分钟作为时间t 的单位。在0

=[（Z/120000）（2000/）*+C] =Z/432+C* 又因为：F(0)=0 所以：C=-Z/432 所以：y=Z/432[1- ] 求得以特解为： F（t）= Z/432[1- ] 在0

数学建模实验报告

数学建模实验报告

一、实验目的 1、通过具体的题目实例，使学生理解数学建模的基本思想和方法，掌握 数学建模分析和解决的基本过程。 2、培养学生主动探索、努力进取的的学风，增强学生的应用意识和创新 能力，为今后从事科研工作打下初步的基础。 二、实验题目 （一）题目一 1、题目：电梯问题有r个人在一楼进入电梯，楼上有n层。设每个 乘客在任何一层楼出电梯的概率相同，试建立一个概率模型，求直 到电梯中的乘客下完时，电梯需停次数的数学期望。 2、问题分析 （1）由于每位乘客在任何一层楼出电梯的概率相同，且各种可能的情况众多且复杂，难于推导。所以选择采用计算机模拟的 方法，求得近似结果。 （2）通过增加试验次数，使近似解越来越接近真实情况。 3、模型建立 建立一个n*r的二维随机矩阵，该矩阵每列元素中只有一个为1，其余都为0，这代表每个乘客在对应的楼层下电梯（因为每 个乘客只会在某一层下，故没列只有一个1）。而每行中1的个数 代表在该楼层下的乘客的人数。 再建立一个有n个元素的一位数组，数组中只有0和1,其中1代表该层有人下，0代表该层没人下。 例如： 给定n=8;r=6（楼8层，乘了6个人）,则建立的二维随机矩阵及与之相关的应建立的一维数组为： m = 0 0 1 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0 1 0 0 c = 1 1 0 1 0 1 1 1 4、解决方法（MATLAB程序代码）：

n=10;r=10;d=1000; a=0; for l=1:d m=full(sparse(randint(1,r,[1,n]),1:r,1,n,r)); c=zeros(n,1); for i=1:n for j=1:r if m(i,j)==1 c(j)=1; break; end continue; end end s=0; for x=1:n if c(x)==1 s=s+1; end continue; end a=a+s; end a/d 5、实验结果 ans = 6.5150 那么，当楼高11层，乘坐10人时，电梯需停次数的数学期望为6.5150。 （二）题目二 1、问题：某厂生产甲乙两种口味的饮料,每百箱甲饮料需用原料6 千克,工人10名,可获利10万元;每百箱乙饮料需用原料5千 克,工人20名,可获利9万元.今工厂共有原料60千克,工人 150名,又由于其他条件所限甲饮料产量不超过8百箱.问如何 安排生产计划,即两种饮料各生产多少使获利最大.进一步讨 论: 1)若投资0.8万元可增加原料1千克,问应否作这项投资. 2)若每百箱甲饮料获利可增加1万元,问应否改变生产计划. 2、问题分析 （1）题目中共有3个约束条件，分别来自原料量、工人数与甲饮料产量的限制。 （2）目标函数是求获利最大时的生产分配，应用MATLAB时要转换

数学建模与数学实验试卷及答案

数学建模与数学实验试卷及答案 二、本题10分(写出程序和结果) 蚌埠学院2010—2011学年第二学期 2,x在 [-5 ，5] 区间内的最小值，并作图加以验证。求函数yxe,，,3《数学建模与数学实验》补考试卷答案 f1=inline('x.^2 +exp(-x)-3') 注意事项:1、适用班级:09数学与应用数学本科1,2班 2、本试卷共1页，附答题纸1页。满分100分。 x=fmin(f1,-5,5) 3、考查时间100分钟。 y=f1(x) 4、考查方式:开卷 fplot(f1,[-5,5]) 一、填空:(每空4分，共60分) x = 0.3517，y== -2.1728 123111,，，，， ,,,,三、本题15分(写出程序和结果) 1. 已知，，则A的秩为 3 ，A的特征值为 A,612B,234,,,, ,,,,,215531,，，，，360000xx，,,12,max2.5fxx,，求解:， stxx..250000，,,1212-1.9766 4.4883 + 0.7734i 4.4883 - 0.7734i ，若令 A([1,3],:)= B([2,3],:)，则,x,150001,A(2,:)= 6 1 2 ; 解: xxx，，,22,123,model: 2. 的解为 1.25 ,0.25 0.5 ; xxx，，,521,123max=2.5*x1+x2; ,242xxx，，,123,3*x1+x2<=60000; 装订线内不要答题 2*x1+x2<=50000; 3. 将1234521 分解成质因数乘积的命令为_factor(sym(‘1234521’)),

数学建模实验

数学建模课程实验报告 专题实验7 班级数财系1班学号2011040123 丛文 实验题目常微分方程数值解 实验目的 1．掌握用MATLAB求微分方程初值问题数值解的方法； 2．通过实例学习微分方程模型解决简化的实际问题； 3．了解欧拉方法和龙格库塔方法的基本思想。 实验容 (包括分 析过程、 方法、和 代码，结 果) 1. 用欧拉方法和龙格库塔方法求下列微分方程初值问题的数值 解，画出解的图形，对结果进行分析比较 解；M文件 function f=f(x,y) f=y+2*x; 程序； clc;clear; a=0;b=1; %求解区间 [x1,y_r]=ode45('f',[a b],1); %调用龙格库塔求解函数求解数值 解； %% 以下利用Euler方法求解 y(1)=1;N=100;h=(b-a)/N; x=a:h:b;

for i=1:N y(i+1)=y(i)+h*f(x(i),y(i)); end figure(1) plot(x1,y_r,'r*',x,y,'b+',x,3*exp(x)-2*x-2,'k-');%数值解与真解图 title('数值解与真解图'); legend('RK4','Euler','真解'); xlabel('x');ylabel('y'); figure(2)

plot(x1,abs(y_r-(3*exp(x1)-2*x1-2)),'k-');%龙格库塔方法的误差 title('龙格库塔方法的误差') xlabel('x');ylabel('Error'); figure(3) plot(x,abs(y-(3*exp(x)-2*x-2)),'r-')%Euler方法的误差 title('Euler方法的误差') xlabel('x');ylabel('Error');

大气污染论文数学建模

大气污染评价与预报模型 摘要 本文对空气质量的评价及污染预报问题进行了分析，运用层次分析法依据处理后的数据对六个城市的空气质量进行了具体细致的排序；对2010年9月15日至9月21日的各项污染物浓度、各气象参数运用一元多项式回归模型进行了预测；就气象参数所属城市问题及污染物浓度与其的关系建立了相关性分析模型和多元线性回归模型；最后，根据建模过程和结果，我们对相关部门提出了几个具体的建议。 通过将数据附件所给有效数据，即日污染物浓度，转化为对应的月污染物浓度的均值，根据各城市月均污染浓度做出其随时间的走势折线图，分析了各个城市2SO 、2NO 、PM10之间的特点。我们拟根据API 指数值，以二级达标次数为准，对各城市之间的空气质量进行排名，但由于依据API 的区分空气质量等级时灵敏度较低，故采用了层次分析法对空气质量进行排名。由于我们采用了全部数据进行排名，而E 、F 数据较少，故只对ABCD 进行了排名。依据层次分析法得出的排名为：A 、B 、D 、C 。 为了精确预测各城市短期内的数据，本文选用一元多项式回归模型。对2010年的数据进行分析整理，依据回归模型得出其与时间的关系，得出预测值，并得出其置信度为95%的置信区间，结果显示模型的预测效果尚能接受，能够对所要预测数据进行预测。但由于F 城市数据缺失，根据假设做了合理的定性分析，并未对其进行定量预测。 分析空气质量与气象参数之间的关系时，首先根据数据完整性，气象参数应只属于其中一个城市，排除了D 、E 、F 的可能性，再根据相关性分析的方法，确定了气象参数属于A 城市。根据污染物与气象参数之间的因果关系，建立了多元线性回归模型，由于季节对污染物的浓度存在影响，分季节得出各污染物与各气象参数之间的相关系数，定性分析该相关系数，得出污染物与气象参数之间的关系。最后对该系数的理论与实际意义做了检验。 根据以上分析及结果，确定部分与空气质量控制相关的部门，针对其职能提出了诚恳建议。 关键词：API 评价模型 层次分析 一元多项式回归模型 相关性分析 多元回归

数学建模与数学实验报告

数学建模与数学实验报告 指导教师__郑克龙___ 成绩____________ 组员1：班级______________ 姓名______________ 学号_____________ 组员2：班级______________ 姓名______________ 学号______________ 实验1.(1)绘制函数cos(tan())y x π=的图像，将其程序及图形粘贴在此。 >> x=-pi:0.01:pi; >> y=cos(tan(pi*x)); >> plot(x,y) -4 -3 -2 -1 1 2 3 4 -1-0.8-0.6-0.4-0.200.20.40.60.8 1 （2）用surf,mesh 命令绘制曲面2 2 2z x y =+，将其程序及图形粘贴在此。（注：图形注意拖放，不要太大）（20分） >> [x,y]=meshgrid([-2:0.1:2]); >> z=2*x.^2+y.^2; >> surf(x,y,z)

-2 2 >> mesh(x,y,z) -2 2 实验2. 1、某校60名学生的一次考试成绩如下:

93 75 83 93 91 85 84 82 77 76 77 95 94 89 91 88 86 83 96 81 79 97 78 75 67 69 68 84 83 81 75 66 85 70 94 84 83 82 80 78 74 73 76 70 86 76 90 89 71 66 86 73 80 94 79 78 77 63 53 55 1)计算均值、标准差、极差、偏度、峰度，画出直方图；2)检验分布的正态性；3)若检验符合正态分布，估计正态分布的参数并检验参数. （20分） 1） >> a=[93 75 83 93 91 85 84 82 77 76 77 95 94 89 91 88 86 83 96 81 79 97 78 75 67 69 68 84 83 81 75 66 85 70 94 84 83 82 80 78 74 73 76 70 86 76 90 89 71 66 86 73 80 94 79 78 77 63 53 55]; >> pjz=mean(a) pjz = 80.1000 >> bzhc=std(a) bzhc = 9.7106 >> jc=max(a)-min(a) jc = 44 >> bar(a)

重金属污染源的数学建模

承诺书 我们仔细阅读了中国大学生数学建模竞赛的竞赛规则. 我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。 我们知道，抄袭别人的成果是违反竞赛规则的, 如果引用别人的成果或其他公开的资料（包括网上查到的资料），必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出。 我们郑重承诺，严格遵守竞赛规则，以保证竞赛的公正、公平性。如有违反竞赛规则的行为，我们将受到严肃处理。 我们参赛选择的题号是（从A/B/C/D中选择一项填写）： A 我们的参赛报名号为（如果赛区设置报名号的话）： 所属学校（请填写完整的全名）：中国人民解放军国防科学技术大学参赛队员(打印并签名) ：1. 李腾骥 2. 袁雪强 3. 赵泉朴 指导教师或指导教师组负责人(打印并签名)：杨文强 日期：2011 年9 月9 日赛区评阅编号（由赛区组委会评阅前进行编号）：

编号专用页 赛区评阅编号（由赛区组委会评阅前进行编号）： 全国统一编号（由赛区组委会送交全国前编号）：全国评阅编号（由全国组委会评阅前进行编号）：

基于MCMC的城市表层土壤重金属污染分析 摘要：本文针对某城市城区表层土壤重金属污染分析的问题，从各分区重金属元素污染的空间分布、污染程度、污染原因及传播途径等多方面进行了全面的统计与研究，同时采用贝叶斯推理和马尔科夫链蒙特卡洛法（MCMC），利用MATLAB进行编程，对所设计的模型进行求解，得到了各重金属污染源的地理位置坐标，最后对城市地质环境演变模式的问题进行了探究，得出了初步的结论。 首先，根据所给数据信息，利用MATLAB及EXCEL作图，得到了该地区的地势图及各重金属污染程度的分布图。 对于问题1，通过统计各重金属的重污染点在各区域的分布和计算不同区域各中金属污染程度之间的相关系数，得到了各重金属的大致空间分布，并采用改进的内梅罗指数法，得出不同区域重金属的污染程度。 对于问题2，在问题1所得数据及结论的基础上，结合实际情况，对各重金属污染的主要原因进行了较为全面的分析。 对于问题3，通过对污染的传播特征的研究，基于城市的地质特性，得到空气传播为污染主要传播方式的结论，并由此建立关于气体传播的概率模型，采用贝叶斯推理和马尔科夫链蒙特卡洛法，利用MATLAB进行50000次抽样，得到收敛的坐标值，从而计算出各重金属污染源的地理位置坐标分别为：As （2742,7294）、Cd（2382,3693）、Cr（3470,2309）、Cu（2707,2295）、Hg(2708,2294)、Ni（1647,2729）、Pb（2882,3617）、Zn（4152,2299）。 对于问题4，以题中城市土壤重金属污染的演变规律为例进行了一定深度的探究，并对所需信息及模型设计方案进行了初步的讨论。 关键词：内梅罗指数法相关系数贝叶斯推理MCMC法MH算法 一、问题的提出与分析 （一）背景介绍 随着城市经济的快速发展和城市人口的不断增加，人类活动对城市环境质量的影响日显突出。对城市土壤地质环境异常的查证，以及如何应用查证获得的海量数据资料开展城市环境质量评价，研究人类活动影响下城市地质环境的演变模式，日益成为人们关注的焦点。 按照功能划分，城区一般可分为生活区、工业区、山区、主干道路区及公园绿地区等，分别记为1类区、2类区、……、5类区，不同的区域环境受人类活动影响的程度不同。 现对某城市城区土壤地质环境进行调查。为此，将所考察的城区划分为间距1公里左右的网格子区域，按照每平方公里1个采样点对表层土（0~10 厘米深度）进行取样、编号，并用GPS记录采样点的位置。应用专门仪器测试分析，获得了每个样本所含的多种化学元素的浓度数据。另一方面，按照2公里的间距在那些远离人群及工业活动的自然区取样，将其作为该城区表层土壤中元素的背景值。

《数学建模与数学实验》课程论文

10级信息《数学建模与数学实验(实践)》任务书 一、设计目的 通过《数学建模与数学实验(实践)》实践环节，掌握本门课程的众多数学建模方法和原理，并通过编写C语言或matlab程序，掌握各种基本算法在计算机中的具体表达方法，并逐一了解它们的优劣、稳定性以及收敛性。在熟练掌握C 语言或matlab语言编程的基础上，编写算法和稳定性均佳、通用性强、可读性好，输入输出方便的程序，以解决实际中的一些科学计算问题。 二、设计教学内容 1线性规划（掌握线性规划的模型、算法以及Matlab 实现）。整数线性规划（掌握整数线性规划形式和解法）。 2微分方程建模（掌握根据规律建立微分方程模型及解法；微分方程模型的Matlab 实现）。 3最短路问题（掌握最短路问题及算法，了解利用最短路问题解决实际问题）。 行遍性问题（了解行遍性问题，掌握其TSP算法）。 4回归分析（掌握一元线性回归和多元线性回归，掌握回归的Matlab实现）。 5计算机模拟（掌握Monte-carlo方法、了解随机数的产生；能够用Monte-carlo 解决实际问题）。 6插值与拟合（了解数据拟合基本原理，掌握用利用Matlab工具箱解决曲线拟合问题）。 三、设计时间 2012—2013学年第1学期：第16周共计一周 目录 一、10级信息《数学建模与数学实验(实践)》任务书 (1) 二、饭店餐桌的布局问题 (3) 摘要 (3)

问题重述 (3) 模型假设 (3) 模型分析 (4) 模型的建立和求解 (4) 模型推广 (9) 参考文献 (9) 三、白酒配比销售问题 (10) 摘要 (10) 问题重述 (11) 问题分析 (12) 模型假设 (12) 符号及变量说明 (12) 模型的建立与求解 (13) 模型的检验 (18) 模型的评价与推广 (19) 附录 (21) 饭店餐桌的布局问题 摘要 饭店餐桌的布局对于一个饭店有着很重要的作用。本文讨论的就是饭店餐桌的布局问题，根据实际需求及规定建立模型，同时考虑餐桌的类型及规格，尤其是餐桌的摆放技巧，保证使饭店能容纳的人数达到最大。根据所需餐桌的数量

污染物扩散模型-深圳数学建模

赛区评阅编号（由赛区组委会填写）： 2015高教社杯全国大学生数学建模竞赛 承诺书 我们仔细阅读了《全国大学生数学建模竞赛章程》和《全国大学生数学建模竞赛参赛规则》（以下简称为“竞赛章程和参赛规则”，可从全国大学生数学建模竞赛网站下载）。 我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。 我们知道，抄袭别人的成果是违反竞赛章程和参赛规则的，如果引用别人的成果或其他公开的资料（包括网上查到的资料），必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出。 我们郑重承诺，严格遵守竞赛章程和参赛规则，以保证竞赛的公正、公平性。如有违反竞赛章程和参赛规则的行为，我们将受到严肃处理。 我们授权全国大学生数学建模竞赛组委会，可将我们的论文以任何形式进行公开展示（包括进行网上公示，在书籍、期刊和其他媒体进行正式或非正式发表等）。 我们参赛选择的题号（从A/B/C/D中选择一项填写）： C 我们的报名参赛队号（12位数字全国统一编号）： 参赛学校（完整的学校全称，不含院系名）：温州医科大学 参赛队员 (打印并签名) ：1. 章成俊 2. 杨超 3. 谢锦 指导教师或指导教师组负责人 (打印并签名)： 日期：年月日

赛区评阅编号（由赛区组委会填写）： 2015高教社杯全国大学生数学建模竞赛 编号专用页 送全国评阅统一编号（由赛区组委会填写）： 全国评阅随机编号（由全国组委会填写）：

对垃圾处理厂污染的动态监控及居民补偿 摘要 城市垃圾处理问题是一个世界性难题。目前垃圾焚烧正逐步成为中国垃圾处理的主要手段之一。本论文构根据题目设置的垃圾处理厂规模，建立了环境动态监控体系，并根据潜在污染风险对周围居民进行了合理经济补偿的设计。 对于问题（1），为了实现对垃圾焚烧厂烟气排放及相关环境影响状况的动态监控，本论文在高斯烟羽模型的基础上进行改进，引入温度、降雨对污染物扩散的影响，建立了新的污染物扩散模型。本论文创新性的提出了风雨影响指数M,用来衡量风向、降雨对颗粒物扩散的影响。本论文将抽象的污染物含量形象化，利用空气污染指数API描述具体的污染程度及其给周围居民带来的影响。并且从不同角度给出了模型检验，验证了所建模型的准确性。 对于问题（1）具体赔偿方案的制定，在综合考虑了不同方位风向频率、受污染时间、受污染程度的基础上，本论文使用了层次分析法，并且进行了一致性检验，使得赔偿方案具有说服力。通过MATLAB编程，计算出当政府和垃圾处理厂共支付风险赔偿金为N时，得出居住地的每位居民应得的赔偿金额计算公式。对于监测点的设置，经计算共需21个，具体布置情况见后文。 对于问题（2），在题目所述的发生事故的情况下，对污染物的具体含量进行了合理的预测与假设。模拟出酸性物质与颗粒物的影响范围，并根据具体的污染程度设置不同的污染区。对每个污染区的不同情况设置更改监测点的设置，并且在问题（1）的基础上对居民的经济补偿进行合理修改。 关键词：高斯烟羽模型，层次分析法，空气污染指数，烟气抬升公式 一、问题重述 “垃圾围城”是世界性难题，在今天的中国显得尤为突出。数据显示，目前全国三分之二以上的城市面临“垃圾围城”问题，垃圾堆放累计侵占土地75万亩。因此，垃圾焚烧正逐步成为中国垃圾处理的主要手段之一。然而，由于政府监管不力、投资者目光短浅等多方面的原因，致使前些年各地建设的垃圾焚烧电厂在运营中出现了环境污染问题，给垃圾焚烧技术在我国的推广造成了很大阻力，许多城市的新建垃圾焚烧厂选址都出现因居民反对而难以落地的局面。在垃圾焚烧厂运行监管方面，目前主要是在垃圾焚烧厂内进行测量监控，缺少从周边环境视角出发的外围动态监控，因而难以形成为民众所信服的全方位垃圾焚烧厂环境监控体系。 深圳市某地点计划建立一个中型的垃圾焚烧厂，计划处理垃圾量1950吨/天（设置三台可处理垃圾650吨/天的焚烧炉，排烟口高度80米，每天24小时运转）。从构建环境动态监控体系、并根据潜在污染风险对周围居民进行合理经济补偿的需求出发，有关部门希望能综合考虑垃圾焚烧厂对周围带来环境污染以及其他危害的多种因素（例如，焚烧炉的污染物排放量、居住点离开垃圾焚烧厂的距离、风力和风向及降雨等气象条件、地形地貌以及建筑物的遮挡程度等等），在进行科学定量分析的基础

数学建模与数学实验

数学建模与数学实验 实验报告 班级: 数学师范153 姓名:付爽 学号:1502012060 实验名称: 数列极限与函数极限 基础实验 基础实验一数列极限与函数极限第一部分实验指导书解读

一、实验目的 从刘徽的割圆术、裴波那奇数列研究数列的收敛性并抽象出极限的定义;理解数列收敛的准则;理解函数极限与数列极限的关系。 二、实验使用软件 Mathematic 5、0 三.实验的基本理论即方法 1割圆术 中国古代数学家刘徽在《九章算术注》方田章圆田术中创造了割圆术计算圆周率π。刘徽先注意到圆内接正多边形的面积小于圆面积;其次,当将边数屡次加倍时,正多边形的面积增大,边数愈大则正多边形面积愈近于圆的面积。 “割之弥细,所失弥少。割之又割以至不可割,则与圆合体而无所失矣。”这几句话明确地表明了刘徽的极限思想。 以n S 表示单位圆的圆内接正1 23-?n 多边形面积,则其极限为 圆周率π。用下列Mathematica 程序可以从量与形两个角度考察数列{n S }的收敛情况: m=2;n=15;k=10; For[i=2,i<=n,i++, l[i_]:=N[2*Sin[Pi/(3*2^i)],k]; (圆

内接正1 23-?n 多边形边长) s[i_]:=N[3*2^(i-1)*l[i]*Sqrt[1-(l[i])^2/4],k]; (圆内接正1 23-?n 多边形面积) r[i_]:=Pi-s[i]; d[i_]:=s[i]-s[i-1]; Print[i," ",r[i]," ",l[i]," ",s[i]," ",d[i]] ] t=Table[{i,s[i]},{i,m,n}] (数组) ListPlot[t] (散点图) 2裴波那奇数列与黄金分割 由2110;1; 0--+===n n n F F F F F 有著名的裴波那奇数列}{n F 。 如果令n n n F F R 11 --=,由n F 递推公式可得出 11111/11---+=+=+=n n n n n n n R F F F F F R ,]251251[511 1 ++??? ? ??--??? ? ??+=n n n F ; 2 15lim lim 1 -==+∞ →∞ →n n n n n F F R 。 用下列Mathematica 程序可以从量与形两个角度考察数列{n R }的收敛情况: n=14,k=10; For[i=3,i<=n,i++, t1=(Sqrt[5]+1)/2; t2=(1-Sqrt[5])/2;