基于差分方程的湖水污染问题

基于差分方程模型对湖水污染问题的分析一 概要 本文讨论了湖水浓度变化趋势的预测问题。 通过分析水流输入输出湖泊的过程，建立了湖水污染浓度随时间变化的含参变量的微分方程模型，在湖水污染浓度恒定和自然净化速率呈线性关系的情况下，求得其准确解，带入具体数据得到结论：在PCA声称的湖水污染浓度下，湖的环境不会恶化；在工作人员实施测得的河水浓度下，湖的环境将会恶化。同时建立了计算机模拟模型，带入具体数值，运用时间步长法来仿真模拟了在湖水污染浓度稳定以前湖水每天的变化情况，输出自PCA建厂以来每年的湖水污染浓度变化，得到与微分方程模型相同的结论。在全停产和半停产时，通过前面的两个模型可以计算湖水污染浓度在自然进化影响下的恢复到自然净华指标所需年限。并可得到结论：在半停产条件下，在选定的自然净化速率常数的约束下，只有当河水污染浓度降至原来的3.15%（自然净化速率成线性关系），4.7%（自然净化速率呈指数关系），才有可能使河水在100年内恢复至0.001，然后给出整改建议。关键字：时间步长法 微分方程模型 logistic模型二 问题重述Pure河是流入Pristine湖的唯一河流。50年前PCA公司在此河旁建起一个生产设施并投入运行。PCA将为处理的废水排入河中,导致了Pristine湖被污染。PCA公司声称:已排放的废水的标准多年从未改变且不会对湖的环境有影响。Pristine湖的湖容量为，流入(流出)的水流速度为年。PCA声称河水污染浓度仅为0.001mol/l，自工厂开工以来没有改变过。问题：1.在花费时间和经费去测试之前，建立数学模型用PCA提供的公开数据判断湖的环境是否会恶化。2.派出野外工作人员测得目前湖水污染浓度为0.03 mol/l，再测得河水污染浓度为0.05 mol/l。以新数据为依据考虑湖水污染问题的数学模型。3.现在假设你是环保局的所聘请的高级顾问，请向你的雇主提交一份报告.内容包括:(1)在工厂停产(或半停产)条件下，湖水自然净化所需年限(净化指标为污染浓度不超过0.001 mol/l)；(2)为保护环境，对PCA进行整改的建议。三 问题假设1.假设降水量和蒸发量相等。2.湖中流入量和流出量相等且一直未变。3.湖水混合均匀，污染物质和湖水瞬间混合。4.湖内无其它污染源，只有PAC的污染源。5.在整个大环境下，河水和湖水在污染源的浓度不一致。6.在时间允许的范围内，假设整个湖水容量不发生变化。7.污染源属于什么态的物质四 符号说明和名词解释 河水污染浓度湖水污染物浓度V湖泊容量Lc自然净化速率流入（流出）的水流速度t从PCA建厂至考察时刻的时间段天五 问题分析5.1 对于水域基本情况分析 由于题目假设可得，我们对于整个水质的大环境假设为一个容积不变的大容器。基于在该流域里，降水量和蒸发量保持相等，湖水流入量和流出量相等且一直不变，所以我们假设该水域环境下，其容积不发生变化。针对污染源的不同污染和净化的情况，我们把该水域分为两个部分，湖水和河水。这两部分水的成分只在湖水的连接区域保持联系。鉴于河水的流动性，水质轮新的速度较快，所以我们可以做这样的定义：其中湖水作为净化的场所，而对于河水而言，其功能上就是作为PC A公司投放污染源的场所。所以我们考虑对于污染源进入该流域内的速度，根据题意直接所得的相关表示形式。5.2 对于题设条件下对于水质净化作用的分析而对于该水域内自带的水质净化作用，我们初步定义，对于水质的净化作用需要考虑该水域的面积和深度（即容积的大小）等因素。对于该题而言，该流域的净化作用能力主要体现在流出过程中的排污能力和多种形式的自我净化作用能力。针对水质的自我净化作用和污染作用这两个过程，我们考虑是一个相互抗衡的阶段，根据题设已知污染源进入水质的速度恒定，我们考虑得，水质的净化作用速度收到该水域环境下的污染浓度限制，这样的过程是一个负反馈的过程。由基本假设得到，只有当达到污染源的浓度超过某个阶段，也就意味着水质的净化作用相应有一定幅度地减少。而针对在无污染的情况下，水质的净化作用随着水质的逐渐改善不断提高。5.3 对于净化标准划定的分析根据题意已知，当该水域的污染程度达到0.001mol/l，才达到需要治理水质的要求，这是一个基本红线。5.4 对于污染源划定的基本分析污染源的划定，我们假设其作为含有某系列物质的溶液，其主要分为含有不同污染能力的固态物质若等，对于题三的整改部分，我们可以对公司提供一定的建议。鉴于公司不能直接切实有效对污染源的量做到减少，只有在排放的部分中做一定改善。对于不同污染能力的污染源物质，可以采取固态物质粉碎颗粒化，对于溶液进一步稀释。综上所述，对于该题，我们得到如下基本分析。在前假设条件的基础之上，湖水容量不变，出河水外无其他的污染源，故我们可将湖泊作为一个封闭的生态系统，其简化的湖水被污染的动态过程为：受污河水输入湖水中，河水与湖水均匀混合，受污水质在湖水区域内进行自我净化，湖水从输出端流出湖泊。湖中污染物的量直接决定了湖水污染浓度，而污染物的量受到以下两方面因素的影响： 1.河水的污染浓度及其流入速度（根据已知此速度不变）， 2.湖水的自我净化速度，前者使其增加，后者使其减少。问题一、二、三的实质都是要分析污染浓度的变化趋势，而关键就在于前一因素的不同。问题一中，河水污染浓度不变，恒为0.001mol/L；问题二中，河水污染浓度可能会变化，受PCA效益的影响而按一定规律波动；问题三中，在全停产或半停产的情况下，和硕污染浓度为0或减为问题二中的一部分。后一因素（自然净化速度）在三个问题中的作用都是相同的。六 模型建立及求解 问题一：我们对于整个水域的探究，假设在整个大环境下，对于该水域而言，基于河水的流动性的特点，我们可以做基本假设，认为PC A公司将污染源投放到河水中，也就是说河水具有污染源进入水域的特点。基于湖水的相对闭塞性的特点，我们可以做基本假设，认为湖水区域作为自我净化的场所。我们可以在该区域找部分区域作为其中净化的场所。在对于整个大环境的探究过程中，我们暂定认为，该水域由一条流进湖的pure 河中间的Pristine湖和一条流出湖的组成。对于该水域的动态变化过程可以认为，该水域的污染来源是PCA公司的污染投放，该水域对于污染的净化能力分为其湖水的自我净化能力和流出该湖的清洁能力。根据微积分的知识可知，在适当短的时间段之内，通过建立微分方程，可以将连续的过程离散化，从而可以得到湖水污染浓度与时间之间的关系式。所以鉴于湖水的实际情况，我们可以用如下图的示意图1表示这种情况。 图1由相关研究表示，湖水的自我净化作用相对于污染物直接从输出端流出的量较少，我们在处理问题一时，先不考虑湖水区域的自我净化作用。（1）不考虑湖水的自我净化作用由于假设湖中流入量和流出量相等，而在经过与湖水均匀混合后，流出湖水污染浓度明显减小，故流出污染物的量小于流入污染物的量，污染物将在湖中沉积，从而使湖水污染浓度增加，当其增加至和输入的河水污染浓度相等时。在这种情况下，湖水污染浓度达到最大，并稳定在这一数值，不再增加。建立湖水污染浓度随时间变化的微分方程模型：设在极短时间dt，湖水污染浓度 增加d，在将湖水被污染这一连续动态过程简化为离散的瞬间静止状态（如问题分析中所述）之后，根湖中剩余量=输入量输出量，我们可以列微分方程如下我们做如下规定： （1）把该式子求解得， （2）将相关数据带入化简可得：可得和t的关系式如下： （3）通过此关系式我们可知，当t=时，湖水污染浓度将趋近与0.001mol/L，即湖的环境不会恶化。 建立计算机模拟模型：湖水污染浓度的变化时有湖中污染物随时间的积累而引起的，这个逐步积累的过程我们可以用计算机进行仿真模拟，其实质为完成一个循环累加的过程，并可改变时间步长，如一年一年的累积，一月一月的累积，一天一天的累积，从而使我们的模拟值逐步精确，可与微分方程求得精确解比较，分析误差。为提高模拟结果的精确性和运算的效率，我们采用了逐天累加，数出年污染浓度的方式。模拟程序见附件一，从后面的分析中我们可知：在此河水污染浓度恒定，无自然净化的最简单的情况下建立的模型是以后问题的基础，后面的问题只是改变条件或数据，其实质是不变的，故我们在此程序中加入了多个选择语句，在不同的条件或数据下执行不同的命令，从而用一个程序解决了全部的模拟问题。模拟所得数据如下表1所示：年份135810 1.734.356.147.828.51年份1520 304050 9.439.789.979.9959.999 表1根据模拟数据所作的湖水污染浓度变化趋势如下图2所示 图2（1） 考虑湖水的自我净化作用自我净化作用的原因：根据资料显示：湖水中的污染物可分为有机污染物和无机污染物两大类，在多种环境因素（阳光，空气，水，水中生物，水中化学物质，重力等）的作用下，通过物理沉降，化学反应和生物转化一系列复杂的活动，它们的量会发生变化。有机污染物在水环境中的迁移过程过程如下图3，4所示。 图3图4可以看出整个过程是相当的复杂，不仅过程多，而且在相同的过程中，不同的物质有着不同的结果。此问题中没有明确给出输出湖中的污染物的种类，也没有对湖泊环境做任何描述，无疑给问题的解决增加了极大的难度。为使问题简化，我们从一般情况出发，假设：湖中所进行的反应均为一级反应，有机物的存在不会对环境参数造成改变，环境固定（自然净化速率恒定），考虑湖泊生态系统中起主要作用的几种过程作简要分析。1.物理沉降。不同物质有着不同的沉降速度常数，沉降速率为2.挥发。不同物质有着不同的挥发速率常数，在有机物在水体上的大气中的分压为0的条件下，挥发速率为 (Z为水体深度) 3.水解反应。不同物质有着不同的水解速率常数，在一级反应的条件下，水解速率为； 4.生物降解反应。不同物质有着不同的降解速率常数，在一级反应的条件下，生物降解速率；观察以上各式可发现：在一定条件下，各个速率与湖中污染程度呈线性关系，又因为某种化学物的总消失速率等于各过程消失速率之和，故我们可以认为：自然净化速率（为正常数，与环境及污染种类有关）如上所述，我们可得上述表达式的改进。 （3）把该式求解得 （4）当时间t=时，即湖水污染浓度稳定于一个与K有关的值。B.在最简模拟程序的基础上，从每天的积累量中减去每天的自然净化量Kx时间步长，重复循环可得到逐年湖水污染度的值，设K=0.0003/天。K越大，湖水污染将越快地稳定于一个更小的浓度值。当=0.001时，若要稳定于0.001，则K值为0，而根据实际情况，任何生态系统都存在一定的自然净化能力，故K应是一个大于0的正实数，将稳定于一个小于的值。这说明，在PC A声称的条件下，湖水污染的浓度将越快地稳定于一个更小的浓度值。根据所得数据如下表2所示年份135810 1.5873.4374.3154.8784.984年份1520304050 5.0875.1035.1065.1065.106表2根据模拟数据所作的湖水污染浓度变化趋势如图5所示， 图5C 进一步考虑自然净化速率的影响。更贴切于实际的情况是，自然净化速率c与湖水污染程度成指数关系，c随的增加而增加，但增加的速率会逐步减小，用关系式表达即可为：， 为环境和污染程度种类有关的常数。将此关系式代入微分方程（3），得到一个新的微分方程，此方程无解。但是我们可以通过此关系式，建立计算机模拟求得数值解。的大小影响着达到最终稳定浓度，的大小影响着达到最终稳定浓度的快慢。在前面达到线性关系的基础上，在=0.001的c值大小确定A至数量级，并稳定年限在10-20年间，从而得到一股估计值。=，=400模拟所得数据如下表3所示年份1358101.4932.9713.5393.8163.87115203040503.9033.9063.9063.9063.906表3根据模拟所得数据所做的湖水污染浓度变化趋势如下图6所示： 图6观察数据可知：无论自然净化速率c与湖水污染浓度呈线性关系还是指数关系，当湖水污染浓度=0.001时，都将稳定于一个小于0.001的值，也即：湖的环境不会恶化。问题二：根据实际情况，我们在此只考虑存在自然净化的情况。（1）河水污染浓度恒定为0.05mol/l，沿用问题一（2）的微分方程模型和差分模拟模型，我们可以得到以下结果：A. 微分方程模型（自然净化速率c与湖水污染浓度呈线性关系）：当时间t=时，即湖水污染浓度稳定于一个与K有关的值。当=0.05mol/l时，若要稳定于0.001mol/l，则K的值为。也即是说:只有当自然净化速率常数K， 湖的环境才不会恶化，最终稳定在小于或等于0.001mol/l的某一水平。根据我们的设定的K=0.0003（）我们可计算，湖水的最终污染浓度=0.0317mol/l,超过净化指标0.001mol/l，故在此条件下，湖的环境将会恶化。A. 差分模拟模型：在c与成线性关系时，得到模拟数据如下表4所示:年份1358100.821.882.462.883.01年133.163.173.173.17 表4在c与成 指数关系时，得到模拟数据如下表5所示：年份1358100.812.012.823.594.02年324.484.564.584.58表5观察以上数据可知：由于河水污染浓度过大，湖水污染浓度从第一年就超过了净化指标，并逐年增加，使湖的环境恶化。（2）.河水污染浓度变化：根据实际情况，一个工厂的生产量并非恒定不变，每年每月甚至每天也有所不同，从而起其排放污染量也将随生产量的变化而变化的。假设PCA自建立以来的年排污染量服从logistic模型，并非考虑到污染物的量受到多方面因素的影响，在每一个时刻的量上加上一个服从正态分布，范围在此量的10%的较小量，则可以建立河水污染浓度的模型如下：（M ,c ,r为与工厂效益有关的常数）设此时工厂处于稳步上升的发展阶段，其变化曲线如下图7所示： 图7通过计算机模拟产生随机数（如图中星点所示），带入模拟程序。在c与成线性关系时，得到的模拟数据如下表6所示：年份1510202235400.1850.8011.7927.2729.58153.0497.13年份508290100112140150260.7625.0633.9617.2626.3628.5632.5表6根据模拟所得数据所做的湖水污染浓度变化趋势图如下图8所示： 图8在c与成指数关系时，得到的模拟数据如下表7所示：年份1510202435450.1690.5861.2505.1949.04644.2695.75年份5080100120130140150321.7925.5959.6950.6949.6952.5964.8表7根据模拟所得数据所作的湖水污染浓度变化趋势如下图9所示：图9 在第50年时，无论自然净化速率c与湖水污染浓度成线性关系还是指数关系，可以看到湖水污染浓度都接近工作人员实地测得值，这从检验的角度说明我们的对模型参数的估计是可取的。观察模型模拟所得数据可知，湖的环境将会恶化。问题三：（1）.在工厂全停产的条件下，河水污染浓度=0mol/l，湖水污染浓度将会在自然净化的作用下由初值0.03mol/l逐年降低。A. 当自然净化速率c与湖水污染浓度呈线性关系时，其微分方程模型为：，湖水自然净化至0.001mol/l，所需年限为t=11.36年。B. 利用计算机循环模拟，其本质与前面的差分模型相同，只是前面是增加，而这里是负增加（污染物的量逐渐减少），同样可改变时间步长使结果趋于准确，我们仍采用逐天累加，输出年湖水污染浓度的模式，当湖水污染浓度减至0.001mol/l时，停止累加，输出累加年限，此年限即为湖水自然净化所需年限。在c与成线性关系时，得到的年限t=12年119天在c与成指数关系时，得到的年限t=12年8天（2）.在半停产状态下，当湖水污染浓度恒为0.005mol/l,设定湖水污染浓度降至净化指标的最大年限为100年，利用“二分法”编程搜索能够保证湖水在100年内达到净化标准的河水污染浓度，所得数据应满足两个条件：1.在100年内达到某一稳定浓度值，2.此浓度值不超过净化指标0.001mol/l根据程序得，我们可以得到，若要使湖水在100年内恢复至0.001mol/l，在自然净化速率成线性关系的条件下，河水的最大污染浓度为原浓度的3.15%，恢复至净化指标的年限t=36年321天，在自然净化速率成指数关系的条件下，河水的最大污染浓度为原浓度的4.77%，t=31年256天。（3）整改方案：a. 根据前面的计算可知，当自然净化速率一定时，若要湖水在一定年限内恢复至净化指标，河水污染浓度应控制在一定的范围以内。自工作人员检测之日起开始整改，在自然净化速率成线性关系的条件下，河水污染浓度满足，工厂排污量百分比满足，恢复年限满足（年），在自然净化速率成指数关系的条件下，。设环保局要求在整改之后湖水污染浓度恢复至净化年限为，通过二分法模型，改变最大年限为，可求得河水污染浓度及生产百分比，从保护环境的角度出发，应要求PCA工厂立即将排污量降至原来的。例如，当=20年时，=2.84%，实际恢复年限为19年334天，即：若环保局要求在20年之内使湖水污染浓度恢复至净化指标，则工厂排污量应立即降至原总量的2.84%。在自然净化速率呈指数关系的条件下，同理。b. 工厂可通过降低产量和进行污水处理来降低排污量。根据所采用的技术措施的作用原理和去除对象，废水处理方法可分为物理处理法，化学处理法和生物处理法三大类。物理方法有隔栅与筛网，沉淀法，气浮法，离心分离等。化学方法有中和法，混凝法，化学沉淀法，氧化还原法，吸附法，离子交换法，膜分离等。生物处理法有好氧生物处理法，厌氧生物处理法，自然生物处理法等。具体采用哪种，PCA厂应根据本厂的经济状况，地理环境条件，废水的本身性质等选取经济有效的方法。七 模型的进一步优化7.1 自然净化速率的改进：资料显示，污染物的很多生物降解过程是通过好氧反应实现的，当污染浓度较大时，需氧量增加，而湖中氧气总量不变，自然净化速率会有所减小，所以，自然净化速率的函数可作为分段函数考虑。前一段还是用前面的指数模型，后面一段可以考虑。模型这样改进后，更加符合实际情况。7.2 河水污染浓度变化规律模型的改进：根据价值定律，一个工厂的发展总是有起有落，其生产总量绕某一中心值波动，故产生河水污染浓度变化应该是一个具有一定周期性，同时有所增加。可以考虑在原来logstic模型后加上一个傅立叶函数（其周期就是工厂的生产周期）和一个正态的随机变量。八 模型的优缺点(1) 对于首先建立的含参数的微分方程模型，在其有精确解的时候，能给出湖水污染浓度随时间变化的规律，得到很好的预测效果，但微分方程并不是任何条件下都有解，在本问题中，只有当河水污染浓度恒定且自然净化速率成线性关系时，微分方程可解，故适用的范围受到了限制。（2） 差分模拟法，即用数值方法来求解微分方程，它解决了微分方程无准确解其情况，适用范围很大。使用一个兼顾运算效率和模拟结果精度的时间步长是该方法的关键。在处理本问题时，采用的时间步长为天。这样处理的精度大约是0.0001,与准确解相当接近。（3）在考虑自然净化速率和河水浓度变化速率时，我们引入了参数。自然净化速率：1 ：，K为正常数，与环境即污染物种类有关 2 ：，A，B为与环境和污染物种类有关的常数。参数的选择对模型结果有较大的影响，从问题二的第一种情况（河水污染浓度恒定）的计算可看出，K存在着9.31这样一个临界值，当K大于此值时，湖水在现有河水污染浓度下是不会恶化的，若要得到准确结果，就需要根据湖的环境和污染物种类测量K的具体指，在此我们利用了题中所给的实地测量的湖水污染物的种类测量K的具体值，在此我们利用了题中所给的实地测量的湖水污染浓度0.03mol/l对K进行估计。通过模拟数据我们可以看到。第五十年时，湖水污染浓度为0.0317mo/l，与实际情况接近，故可认为我们的参数选取是较为合理的。同理于A，B。河水浓度变化规律：（M ,c ,r为与工厂效益有关的常数）工厂的排污量对湖水污染浓度的变化有着较大的影响，尤其是在工厂在全停产和半停产时，我们可以看到，工厂排污量的较小变化就可引起恢复年限的较大增加，因此，参数的估计是很重要的。我们在设定工厂此时在处于稳步上升阶段的条件下，对参数的估计，而工厂实际的发展情况需要实地调查，估计值不可能完全准确，但通过模拟数据，我们可以看到，在第五十年时，无论自然净化速率c与湖水污染浓度呈线性关系还是指数关系，湖水污染浓度都接近工作人员实地测量值（0.03mol/l），这从检验的角度说明我们的对模型参数的估计是可取的。参考