农作物用水量预测及智能灌溉方法

农作物用水量预测及智能灌溉方法 摘 要 针对问题一，画出题中表1数据的散点图后，发现近似抛物线，为从而对需水量为自变量、产量为因变量分别建立线性、二次和三次函数模型，并通过SPSS软件进行拟合和进一步验证，从而得到二次水分生产函数，表达式为:Y=5.186W2?0.007W?478.472，且拟合残差极小。 针对问题二，首先对相关数据进行 标准 excel标准偏差excel标准偏差函数exl标准差函数国标检验抽样标准表免费下载红头文件格式标准下载 化处理，排除量纲的影响，然后建立基于BP神经网络的作物的水分生产函数模型。通过Matlab得到标准化后的各阶段的蒸发蒸腾量与最终产量之间的关系表达式: ==?5.4743 ???? 其中 net=2.7574?ET+9.2319?ET+4.7184?ET+1.7696?ET??????10.4710。 再将所建立的神经网络模型与附录中两种常见的机理模型进行误差 分析 定性数据统计分析pdf销售业绩分析模板建筑结构震害分析销售进度分析表京东商城竞争战略分析 和对比， 寻求最优的作物水分生产函数模型。 针对问题三我们建立了动态规划模型来求解。结合表2晚稻蒸发蒸腾量与产量数据，在总供水量分别为充分灌溉总水量的30%，70%，90%的情况下，确定了动态规划的阶段变量、决策变量、状态变量、系统方程、目标函数、初始条件及约束条件，并选择了Jensen模型和Blank模型作为不同的目标函数，求出各生育期灌水量的分配情况，发现Jensen模型参数比Blank模型的参数对产量的敏感度更高，分配——————————————————————————————————————————————— 方案 气瓶 现场处置方案 .pdf气瓶 现场处置方案 .doc见习基地管理方案.doc关于群访事件的化解方案建筑工地扬尘治理专项方案下载 更为合理。最后利用已知的充分灌溉条件下各阶段的蒸发蒸腾量进行检验，发现各阶段的灌溉水量基本符合，也就是说明该模型的准确性较高。 针对问题四，我们进一步考虑该地区的降水量，重新构建非充分灌溉 制度 关于办公室下班关闭电源制度矿山事故隐患举报和奖励制度制度下载人事管理制度doc盘点制度下载 的优化模型来求解。通过某地近年的降水资料及田间水量平衡方程可重新确定动态规划的相关概念。其中以作物的相对产量最高，灌溉总用水量最小作为目标函数。利用加权法将多目标动态规划转换为单一目标动态规划问题，目标函数方程为: max ????=????????????(1?????)?????，经过多次试验，确定加权系数c为0.65，从而求出不同年份各生育期灌水量的分配情况。 ????1????2????3????4关键词: BP神经网络;动态规划;智能灌溉 1 1 问题重述 1.1 问题的背景 中国的水资源分布不均，矛盾突出。在农业现代化方面，随着水资源供需矛盾的日益加剧，可采用智能灌溉，以农作物实际需水量为依据，提高灌溉精确度，实施合理的灌溉方法，进而提高水的利用率。 1.2 问题的相关信息 灌溉水利用系数是指在一次灌水期间被农作物利用的净水量与水源渠首处总引进水量的比值，是评价农业水资源利用，指导节水灌溉和大中型灌区续建配套及节水改造健康发展的重要参考。 ——————————————————————————————————————————————— 水分生产函数是指作物产量与水分之间存在的一种投入与产出的数学关系。作物水分生产函数的单因子模型中自变量的形式可以为灌水量、实际腾发量、土壤含水量等，因变量的形式可以为作物产量，平均产量，边际产量等。 非充分灌溉是指在灌溉水不能完全满足作物的生长发育全过程需水量的情况下，以作物水分生产函数为理论依据，将有限的水科学合理(非足额)安排在对产量影响比较大、并能产生较高经济价值的需水临界期供水，从而建立合理的水量与产量关系及分配模式，在水分利用效率、产量、经济效益三方面寻求有效均衡，实现经济效益最大化。 附录A为某省水利用系数的历史数据。 附录B为作物生育阶段水分的两种常见的数学模型:乘法模型和加法模型 附录C为非充分灌溉优化模型背景材料及提示，包括了决策变量及约束、水量平衡方程。 1.3 需解决的问题 根据题中所给信息，本文将问题细化为以下四个问题，并建立数学模型进行分析和研究。 问题一:基于题目所给表格中作物1实测需水量与产量的对应数据，建立该作物全生育期的水分生产函数的模型，即总产量与需水量之间的解析关系，给出详细过程及拟合效果。 问题二:作物的全生育期可以分为若干个生育阶段，以水稻为例，可以分为返青、分蘖、拔节孕穗、抽穗开花、乳熟、黄熟6个生育阶——————————————————————————————————————————————— 段。不同阶段灌溉水量不足均会对最终的产量有影响，基于体重所给某地晚稻分蘖至乳熟各阶段受旱情况对产量影响的数据，利用附录中材料，选取某类优化算法，寻求最优的作物水分生产函数模型，得到各阶段的蒸发蒸腾量(可以理解为灌水量)与最终产量之间的关系。给出详细过程，并将所得结果与附录中常见的机理模型作对比。 问题三:基于已有材料对于非充分灌溉的描述，结合表2晚稻蒸发蒸腾量与产量数据，选择合适的作物水分生产函数(可以选择已得到的模型或者机理模型)，建立该地区水稻的非充分灌溉制度的优化模型。在总供水量分别为充分灌溉总水量的30%，70%，90%的情况下，解决供水量在各生育期(从分蘖到乳熟阶段)合理分配的问题，应包含阶段变量、决策变量、状态变量、系统方程、目标函数、初始条件及约束条件等描述(背景材料及提示见附件)，最后给出求解方案。 问题四:题目中的第三张表格列出了某地近年的降水资料。若进一步考虑该地区的降水量，重新构建非充分灌溉制度的优化模型，从而建立在各生育期合理分配水量方法，使得最大化作物产量的同时，灌溉总用水量最小。 2 2 模型假设与符号说明 2.1 模型的假设 假设1:作物的蒸发蒸腾量等于灌水量; 假设2:作物的总产量只和灌水量有关，不受其他因素影响; 假设3:土壤渗透量忽略不计; ——————————————————————————————————————————————— 2.2 模型的符号说明 符号 符号说明 灌水量 W ET???? Y θ 实际腾发量 土壤含水量 作物产量 注:其它符号在文中具体说明。 3 问题分析 3.1 问题一的分析 对于问题一，首先结合表1中作物1的实测需水量与产量的对应数据，画出散点图，通过图像找出需水量与产量的一些相关规律，发现近似抛物线，从而对需水量为自变量、产量为因变量分别建立线性、二次和三次函数模型，并用题目中数据通过SPSS软件进行拟合和进一步验证，从而得到一个最合适的水分生产函数。 3.2 问题二的分析 对于问题二，首先对相关数据进行标准化处理，排除量纲的影响。然后建立了基于BP神经网络的作物的水分生产函数模型。通过Matlab得到标准化后的各阶段的蒸发蒸腾量与最终产量之间的关系表达式。再将所建立的神经网络模型与附录中两种常见的机理模型进行误差分析和对比，寻求最优的作物水分生产函数模型。 3.3问题三的分析 对于问题三，我们建立了动态规划模型来求解。结合表2晚稻蒸发蒸腾量与产量数据，在总供水量分别为充分灌溉总水量的30%，70%，90%的情况下，确定了动态规划的阶段变量、决策变量、状态变量、——————————————————————————————————————————————— 系统方程、目标函数、初始条件及约束条件，其中，选择了已知两种作物水分生产函数作为不同的目标函数，从而求出不同目标函数下各生育期灌水量的分配情况，选取相应的敏感度最高的分配方法作为最优的灌溉制度。最后利用已知的充分灌溉条件下各阶段的蒸发蒸腾量来对两种分配方法进行检验。 3.4问题四的分析 对于问题四，要进一步考虑该地区的降水量，重新构建非充分灌溉制度的优化模型来求解。通过某地近年的降水资料及田间水量平衡方程可以重新确定动态规划的阶段变量、决策变量、状态变量、系统方程、目标函数、初始条件及约束条件。其中以作物的相对产量最高，灌溉总用水量最小作为目标函数。利用加权 3 法将多目标动态规划转换为单一目标动态规划问题，从而求出不同年份各生育期灌水量的分配情况。 4 问题一的解答 4.1 模型的准备 表4.1 全生命周期实测需水量与产量 处理号 需水量(m3/亩) 产量(kg/亩) 1 187.14 247.63 2 238.38 361.18 3 283.12 430.37 4 315.37 462.89 5 337.65 476.88 ——————————————————————————————————————————————— 6 356.21 483.25 7 387.68 483.05 8 414.46 472.01 9 435.57 456.23 10 456.82 434.08 表4.1为作物1实测需水量与产量的对应数据，根据表中的10组数据，利用SPSS软件，以需水量为横坐标，产量为纵坐标，绘制相应的散点图(如下图)。 图4.1全生命周期实测需水量与产量散点图 4.2 全生育期的水分生产函数模型的建立 由图4.1可以看出，十个点的分布大致符合二次函数图像的特征:当需水量小于350????3?亩时，产量随着需水量的增加而增大;当需水量大于400????3? 亩时， 4 产量随着需水量的增大而减少;当需水量位于350????3?亩到400????3?亩之间时，产量会出现一个峰值。 我们设产量为因变量Y，设需水量为自变量W，利用SPSS软件建立二次曲线模型，同时，为了考虑的更加全面，我们对W，Y也进行了线性拟合和三次曲线拟合，并对这三个拟合结果进行检验和误差分析，最终得出最适合的全生育期的水分生产函数模型。 4.3 模型的求解与分析 ——————————————————————————————————————————————— SPSS软件中的相关操作如下:分析-回归-曲线估计，分别选中线性、二次项和立方模型，点击右上角保存。 运行所得结果如下: 表4.3 模型统计及参数评估 因变量:Y 方程式 R平方 0.558 模型摘要 F df1 df2 显著性 10.112 1 8 0.013 5722062 7 0.000 015.8 5722062 7 0.000 015.8 常数 214.075 参数评估 b1 b2 0.635 b3 线性 二次曲1.000 线模型 三次曲1.000 线模型 -478.472 5.186 -0.007 -478.472 5.186 -0.007 0.000 自变量为W。 由表4.1可以得到线性拟合所得的水分生产函数为 Y=0.635W+214.075;二次曲线模型和三次曲线模型所得结果一样，都为Y=5.186W2?0.007W?478.472。 且在二次和三次曲线模型中，它们的显著性为0，F值很大，决定系数R平方都等于1，表示预测值和观测值完全一致。 利用SPSS的绘图功能将三种估计的曲线和观测值都放入一张图内(图4.4)，观察可得线性曲线和观测值相差较大，而二次曲线和三次曲线几乎完全和观察值重合。 5 图4.3曲线拟合图 通过观察比较，最终得到的水分生产函数是 ——————————————————————————————————————————————— Y=5.186W2?0.007W?478.472。 4.4 误差分析 为了对所建立的模型进行误差分析，我们利用表4.1的十组数据 和SPSS软件中的相关操作，得到三种模型的拟合残差。 表4.4 残差 W Y 线性拟合残差 二次拟合残差 三次拟合残差 238.38 361.18 -4.2625 -0.00365 -0.00365 283.12 430.37 36.5183 0.01208 0.01208 315.37 462.89 48.56005 0.00016 0.00016 337.65 476.88 48.40259 -0.0043 -0.0043 356.21 483.25 42.98728 -0.00132 -0.00132 387.68 483.05 22.80432 -0.00736 -0.00736 414.46 472.01 -5.24056 0.00618 0.00618 435.57 456.23 -34.42509 -0.00401 -0.00401 456.82 434.08 -70.06851 0.00351 0.00351 从上述表格中可以明显发现，二次拟合残差很小，说明实际观察 值与估计值(拟合值)之间的差很小，准确性较高。 5 问题二的解答 5.1 模型的准备 作物水分生产函数是描述作物产量与水分因子之间的数学关系。 6 表5.1 某地晚稻蒸发蒸腾量及产量 1分蘖 ?2拔节孕穗 ?3抽穗开花 ?4乳熟 ?产量 处理编号 ——————————————————————————————————————————————— 处理特征 (mm) (mm) (mm) (mm) (kg/hm2) 0 充分灌溉 148.1 111.8 124.7 89.4 7138.5 1轻旱 1 ?113.2 96.6 92.1 67.2 5757.0 1重旱 2 ?107.6 88.3 84.9 64.2 4576.5 2轻旱 3 ?133.9 91.0 106.9 70.3 6111.0 2重旱 4 ?132.1 77.9 93.9 65.0 4555.5 3轻旱 5 ?128.2 99.4 85.3 78.7 5520.0 3重旱 6 ?129.7 92.5 71.9 69.4 5329.5 4轻旱 7 ?140.5 112.9 108.6 68.6 6345.0 4重旱 8 ?135.3 108.0 101.7 65.0 6040.5 1?2中旱 9 ?110.6 83.3 95.2 72.3 5076.0 2?3中旱 10 ?128.4 90.4 83.4 73.6 5442.0 3?4中旱 11 ?130.1 102.6 94.7 61.4 6130.5 表5.1为某地晚稻 分蘖至乳熟各阶段受旱情况对产量影响的数据。首先，为了避免不同量纲带来的麻烦，我们将表格的数据进行如下处理: 1,2,?,????,????????;作物的潜在产量为????????,??????????????2。通过非充分灌溉实验得到作物N设该作物共有t个生育阶段，由表可 得t=4，各阶段的潜在腾发量为(ET????)????,????= 1,2,?,????,??????????????2。将N组实际腾发量和实际产量的数据按 下述公式作标准化处理 啊。 ????????????=(ET????)?????????(ET????)????,????=1,2,?,????,?? ??=1,2,?,???? ????????=(Y????)?????????????,????=1,2,?,???? 组实际 腾发量(ET????)????????,????=1,2,?,????,????=1,2,?,????,????????和实 ——————————————————————————————————————————————— 际产量(Y????)????,????= 5.2基于BP神经网络的水分生产函数模型的建立 一个3层的BP-人工神经网络模型能够实现任意连续映射，3层BP网络模型，如图所示。 式中，????????????，????????一一处理后的作为输入和输出的训练样本，均为比值，无量纲。 7 图5.2 3层BP神经网络模型 图中的????????????表示输入层第????个神经元与隐含层第????个神经元之间的连接权值，????????????表示隐含层第????个神经元与输出层第????个神经元之间的连接权值。设隐含层第????个神经元的阈值为????????，输出层第????个神经元的阈值为????????，则隐含层第????个神经元的输入为 ????????=?????????????????????????????? ????=1n 式中，????????为输入层第????个神经元的输入值;????为输入层的神经元个数。输出为: b????=?????????????? ????=1,2,?,???? ????式中， ????为隐含层神经元个数，????为激励函数，其形式为 ????(1?1+?????????) 其作用是模拟生物神经元的非线性特性。输出层第????个神经元的输入为 ????????=?????????????????????????????? ????=1输出为: ????????=????(????????) 计算过程中，隐含层神经元的激励函数——————————————————————————————————————————————— 取为S型，二输出层神经元的激励函 数取为非线性。为了提高网络的性能，减少其陷入局部极小值的 可能性，提高收敛速度，通常采用改进的BP算法——动量法来实现， 即 ????????????(t+1)=????????????(????)+?????????????(????+1)+?? ???????????????(????) ?????????????=????(?????????????????) 式中，????为动量因子，????为误差函数。时刻????网络的误差函数 定义为: E(t)=1?2??????????(????)?????????(????)? ????=1????2 式中，????????(????)为输出层第????个神经元在时刻t的实际输 出，????????(????)为该时刻的希望输出，????为输出层的神经元数。 当E(t)?????(????为预先给定的误差)时，网络停止训练，此时的 网络模型就是所需要的。 8 5.3 模型的求解与分析 利用Matlab对所建立的模型进行求解，其中 net=???????????????????????????????????????????????????????????? ???????????????????????; ????1????2?????=??????; 34y=1+?????????????????。 ?????????????????????????????????????????????=6.3511。 ——————————————————————————————————————————————— 将运行结果结合模型准备中标准化处理后的公式，最终得到的各阶段的蒸发蒸腾量(可以理解为灌水量)与最终产量之间的关系式为: y==?5.4743 ???? 其中， net=2.7574????1+9.2319????2+4.7184????3+1.7696????4?10.4710。 5.4 模型的比较和误差分析 由题中已知条件可得两种普遍运用的作物生育阶段水分的数学模型: (1)Jensen模型。 ????????=?? ????????????????=1????求解 得:??????????????????????????????????????????????????????????? ??????????????????=(2.7574 9.2319 4.7184 1.7696); 其中，Y????为实际产量;Y????为供水充足条件下作物最大产量;????????为作物生育期 第i阶段作物缺水敏感系数;????????????为作物实际耗水量;????????????为供水充足条件下作物需水量。 (2)Blank模型。 ????????=??????????? ????????????????=1????其中，????????为不同阶段缺水对产量的敏感系数。 表5.4.1 Jensen模型和Blank模型的敏感系数 模型 ?1分蘖 ?2拔节孕穗 ?3抽穗开花 ?4乳熟 0.2090 0.7025 0.2199 0.1523 Jensen模型λ指数 0.0525 0.5575 0.2909 0.0607 Blank模型A系数 ——————————————————————————————————————————————— 根据表5.5.1中 Jensen模型λ、指数Blank模型A系数以及神经 网络模型求解结果，利用各阶段的实际蒸发蒸腾量对这三种模型进行 误差分析，计算所得的误差结果如下表。 9 表5.4.2 神经网络、Jensen和Blank模型的误差比较 编号 Jensen模型相对误差 Blank模型相对误差 神经网络 1 0.052424135 0.029960265 0.002594 2 0.08010986 0.123228485 0.046493 3 0.077564984 0.067412217 0.021427 4 0.062415559 0.094500888 0.043103 5 0.042334111 0.026204321 0.079795 6 0.027863749 0.03281645 0.014704 7 0.04390797 0.026853212 0.019683 8 0.030922393 0.025651654 0.021224 9 0.018283249 0.020655068 0.006431 10 0.025480277 0.028220418 0.014784 11 0.05154368 0.044767903 0.014645 MAPE 0.0466 0.0473 0.025898 通过比较，可以看出数据拟合精度由高到低的模型依次为BP神 经网络模型、Jensen模型和Blank模型。神经网络拟合的精度最高， 说明通过BP神经网络构建的作物水分生产函数模型最优。 图5.4 神经网络拟合情况 ——————————————————————————————————————————————— 从图5.4中可以看出，神经网络预测值和实际观测值的拟合程度较高，说明该模型的效果较好。 6 问题三的解答 6.1 模型的准备 31分蘖 2拔节孕穗 4乳熟 总水量???产量 处理特征 花 (mm) (mm) (mm) (mm) (kg/hm2) (mm) 充分灌溉 148.1 111.8 124.7 89.4 474 7138.5 蒸发蒸腾量可以理解为灌水量。由表6.1可得，在充分灌溉的条件下，某地晚稻分蘖至乳熟4个阶段所需的总水量为474mm。 因此，题目所要求的在总供水量分别为充分灌溉总水量的30%，70%，90%的情况也就是总供水量分别为142.2mm、331.8mm、426.6mm的情况。 合理的非充分灌溉, 将有限的水量在各生育阶段科学地进行限额分配, 使减产损失最小。整个优化灌溉制度设计的过程, 是一个多阶段决策过程, 可采用动态规划模型求解。 6.2 动态规划优化灌溉模型的建立 动态规划是解决多阶段决策问题的一种方法，在每个阶段中都需作出决策，这种连续不断的作出的多阶段决策，构成一个策略。动态规划特点在于， 10 它可以把一个n维决策问题变换为几个一维最优化问题，即把一个优化问题的求解转化为求多个子问题求解，并在求解过程中保存子问题的最优值，避免了子问题的重复求解，从而降低时间复杂度。 ——————————————————————————————————————————————— 动态规划的基本概念: 1)阶段变量:描述阶段的变量称为阶段变量。本题中可以设上述4个生育阶段为阶段变量k=1,2,…,N,N=4。 2)决策变量:表示当过程处于某一阶段的某个状态时，作出的不同决定，从而确定下一阶段的状态。本题中的各个阶段的实际灌水量m????为决策变量。 3)状态变量:表示每个阶段开始所处的自然状况或客观条件。本题中的状态变量为可用于当前及以后各阶段的灌溉水量 ????k=????????? ????=????????k=1,2.?,N 4)系统方程:也就是状态转移方程，用来确定过程由一个状态到另一个状 态的演变过程。若第k阶段状态变量的值????k、该阶段的决策变量????k一经确定，第k+1阶段的状态变量????????+1的值也就确定，即 ????????+1=?????????????k 5)阶段指标函数:用来衡量所实现过程优劣的一种数量指标为指标函数，具有可分离性并满足递推关系。其中，本题的阶段指标函数为供水量????????分配给 ????第k个阶段获得的效益， 用????????(????????)=(????????????????表示。 6)初始条件:????????=????0。 7)约束条件:首先，所有变量都不能为负数。 ——————————————————————————————————————————————— ?状态约束:每个阶段的灌溉量均不能超过总供水量,即0 ??????????????0。 ?决策约束: A).每个阶段的灌溉水量不能超过该阶段最初所拥有的实际水量， 即0??????????????????; B).4个阶段总的灌水量不能超过非充分灌溉相对应的总供水量， 即??????????????????????0; C).非充分灌溉条件下每个阶段得灌溉量都不能超过充分灌溉条 件时对应的阶段的灌溉水量，即?????????????????????????????????。 8)目标函数: 单位面积实际产量????????与最高产量????????的 比值最大化。 ?利用Jensen模型: ????????F1=max?=max ??? ???????????????????????? ???? ?????????????????利用Blank模型: 6.3 模型的求解与分析 根据所建立的动态规划模型，利用Jensen模型时的基本方程为: ?()????????????????=????????????{????????????????(????????)?( ×????????+1?????????????????)},????=1,2,?,N F2=max?=max?????????? ???????????????????? 11 利用Blank模型时的基本方程为: ?????????(????????)=????????????{????????????????(????????)?( +????????+1?????????????????)},????=1,2,?,N 即由k=N 开始，由后向前依次逆推，逐阶段求出最优决策和过 ——————————————————————————————————————————————— 程的最优值 ,直到最后求出????1?(????1)即得到问题的最优解。 为了使程序运算简便，我们将不同阶段的水量都假设成整数，在此条件下进行动态规划，得到该地区水稻的四个阶段在非充分灌溉条件下的分配制度，通过改变目标函数从而得到不同水分生产函数对应的水量分配。 表6.3 不同模型在不同供水量下的水量分配 实际供水量?分蘖?拔节孕?抽穗开?乳熟实际产量 模型 充分灌溉水量(mm) 穗(mm) 花(mm) (mm) 最大产量Jensen模型 22 75 23 22 0.2824 30% Blank模型 0 112 29 1 0.6268 Jensen模型 79 112 83 58 0.7516 70% Blank模型 6 112 125 89 0.9127 Jensen模型 113 112 119 83 0.9260 90% Blank模型 101 112 125 89 0.9463 根据表6.3的数据，可以发现Jensen模型参数比Blank模型的参数对产量的敏感度更高，Blank模型下在总供水量分别为充分灌溉总水量的30%的情况下所得出的分配方案直接忽略了第一阶段的影响即分配灌水量为零，这与实际情况并不完全符合，查阅文献了解到实际情况下出于保守性考虑，担心突发性的干旱气象因素、人为破坏因素等，一般不会选择第一阶段零灌水。因此Jensen模型下的分配方案更为合理。 6.4 模型的检验 利用表6.1充分灌溉条件下各阶段的水量分配的相关数据对Jensen模型进行检验，如下表所示 ——————————————————————————————————————————————— 表6.4 Jensen模型在充分灌溉时的水量分配 1分蘖2拔节孕穗3抽穗开花4乳熟????模型 (mm) (mm) (mm) (mm) Jensen模148 112 125 89 明该模型的准确性较高。将结果和已知材料中的机理模型敏感系数进行联合分析后，发现得出的结果和原模型基本符合，在非充分灌溉条件下要达到产量最大化，灌水量需主要集中在第二第三敏感系数较大即对产量影响较大的两个阶段,可以基本认同该分配方案。 7 问题四的解答 7.1 模型的建立 (1)阶段变量:根据水稻生长过程，将其全生育期划分为N个不等长的生长阶段，阶段变量为n=1,2,?,???? 。 (2)决策变量:决策变量为各生长阶段的实际灌水量????????及实际蒸发量(ET????)????， ????=1,2,?,???? 12 (3)状态变量:状态变量为各阶段初可用于分配的灌溉水量????????及初始田面蓄水深度?????。 (4)系统方程:系统方程是描述系统在运动中状态转移的方程，由于本系统有2 个状态变量，系统方程也有2 个: 第1是水量分配方程，若对第????个生长阶段采用决策????????时, 可表示为: ????????+1=????????????????? 式中: ????????， ????????+1分别为第????及第????+1阶段初系统可用于分配的水量(换算成单位面积的水深); ???????? 为第???? ——————————————————————————————————————————————— 阶段的灌水量。 第2是田间水量平衡方程, 即: ?????+1= ?????+ ????????+ ?????????(ET????)????? ???????? 式中: ????? , ?????+1 :分别为第????及第????+1阶段初田面水层深度; ???????? 为第i 阶段排水量。 (5)目标函数:采用Jensen提出的在供水不足条件下，水量和农作物实际产量的连乘模型，目标函数为单位面积的实际产量????????与最高产量????????的比值最大,以及灌溉总用水量????最小，即: ????????max F= ??=?? ???????????? 由于多目标规划无法直接求解，故本题采取加权求和法将多目标规划问题转换为单一目标规划问题求解。 加权求和法: ????????max ????=????????(1?????)?????=????? ????(1?????)?? ? ??????? ????????????????=1????=1????????????min ????=? ???????? ? ???=1????????=1式中，c为常数，0?c?1。 (6)约束条件 ?决策约束: 0 ? ????????? ???????? ? ????????=???? ????=1????(????????????????????)????? (????????????)????? (????????????????????)???? ????=1 ,2 ,…,N 式中:Q为全生长期单位面积上可供分配的水量; (????????????????????)????、 ——————————————————————————————————————————————— (????????????????????)????分别 为第????阶段的最大及最小蒸发蒸腾量。 ?田面水层 ?????的约束:(?????????????)????? ??????(?????????????)???? 式中: (?????????????)???? , (?????????????)????分别为第????阶段田面水深的上下限。 (?????????????)????的值可以根据实际灌溉经验确定, (?????????????)????则根据根系层土壤平均含水率下限值(根据经验确定，一般取饱和含水率的80%)及根系层深度确定，一般为负值。 (7)初始条件: ?初始田面水深 ?1= ?0，式中 ?0为插秧后稻田保持的水深，一般取20mm。 13 ?第一阶段可用于分配的水量为水稻全生长期可用于分配的水量，即: ????1=????。 7.2 模型的求解与分析 1990年到2010年的降雨量，因此我们要先预测2006年到2010年的水利用系数。 利用spss专家型建模器进行预测，得到2006年到2010 年的预测值为0.53、0.53、0.53、0.53、0.54。得到的拟合预测图如下: 图7.2 水利用系数拟合预测图 水源渠首处总引进水量: ——————————————————————————————————————————————— ????????= ?????+ ????????+ ???????? 水利用系数是指在一次灌水期间被农作物利用的净水量与水源渠首处总引进水量的比值，不妨记为K，即 (ET????)????K= ????????将上式代入田间水量平衡方程，整理得: 有效降雨量的确定:一次降雨不超过田面允许水深上限的降雨量均为有效降雨量，超过部分为排水量。用公式表示为: 14 ?????+1= (1?K)(?????+ ????????+ ????????)? ???????? ???????? ?????????(?????????????)????? (?????????????)???? (?????????????)?????????? ????????>(?????????????)????? ????? 式中: ???????? 为第????阶段降雨量。 对模型求解得: 年份 (mm) (mm) (mm) (mm) 量 2010 118 138 89 1 346 2009 126 128 96 6 356 2008 139 115 108 15 377 2007 134 119 104 12 369 2006 129 124 100 9 362 2005 134 119 104 12 369 2004 134 119 104 12 369 2003 134 119 104 12 369 2002 129 124 100 9 362 2001 134 119 104 12 369 ——————————————————————————————————————————————— 2000 134 119 104 12 369 1999 134 119 104 12 369 1998 139 115 108 15 377 1997 139 115 108 15 377 1996 139 115 108 15 377 1995 142 111 114 17 384 1994 142 111 114 17 384 1993 142 111 114 17 384 1992 142 111 114 17 384 1991 142 111 114 17 384 ????????=?值 0.2360 0.2529 0.2769 0.2690 0.2610 0.2690 0.2690 0.2690 0.2610 0.2690 0.2690 0.2690 0.2769 0.2769 0.2769 0.2836 0.2836 0.2836 0.2836 0.2836 8 模型的评价、改进及推广 8.1 模型的评价 8.1.1模型的优点 问题三中要求给出4个阶段的灌水量，已知给定灌水量，在一般情况下只要求出其中三个阶段即可得到四个阶段的灌溉方案，因此求解的复杂度为o(????3)。 模型三中所运用了动态规划逆序算法则把整个复杂度降低为o(3????2)，有很明显的简化。 问题四中，由于多目标规划无法直接求解，故本题采取加权求和法将多目标规划问题转换为单一目标规划问题求解，极大的简化了问题。 ——————————————————————————————————————————————— 8.1.2模型的缺点 1、忽略了渗透量的影响，虽简化了模型，然而却降低了模型的精确度。 2、为了便于操作，我们只考虑选取Jensen模型时的情况，没有考虑选取Blank模型时的情况，降低了结果的准确性。 8.2 模型的推广 15 我国水资源分布不均匀，而且水资源严重不足，作为一个农业大国，如何在农业中提高水的利用率，达到节水高产的目的对缓解我国水资源供需矛盾具有战略性意义。作物水分生产函数是描述作物产量与水分因子之间的数学关系。无论对节水灌溉的区域规划和系统评估,或是非充分灌溉的推广应用均具有广泛用途。 9 参考文献 [1] 魏占民,陈亚辛,史海滨,等. BP神经网络的春小麦作物-水分模型的初步研究 [ J].灌溉排水, 2002, 21( 2): 12- 16. [2] 崔远来,李远华, 李新健,王春日. 非充分灌溉条件下稻田优化灌溉制度的研究, 水利学报,1995年第10期. [3] 卢文喜,曲武,贾小丰,王喜华. 吉林省西部人工草地优化灌溉制度的研究[J]. 节水灌溉. 2011(01). [4] 叶芳毅,李忠武,李裕元,曾光明,彭亿,张棋,梁红玲,王亚梅. 水稻生长模型发展及应用研究综述[J]. 安徽农业科学. 2009(01). 16 ——————————————————————————————————————————————— 10 附录 附录1:第二问神经网络模型 t=[113.2 107.6 133.9 132.1 128.2 129.7 140.5 135.3 128.4 130.1 96.6 88.3 91 77.9 99.4 92.5 112.9 108 83.3 102.6 92.1 84.9 106.9 93.9 85.3 71.9 108.6 101.7 83.4 94.7 67.2 64.2 70.3 65 78.7 69.4 68.6 65 72.3 61.4] p=[5757 4576.5 6111 4555.5 5520 5329.5 6345 6040.5 5442 6130.5]; % 148.1 111.8 124.7 89.4 7138.5 P=p/7138.5; T(1,:)=t(1,:)/141.8; T(2,:)=t(2,:)/111.8; T(3,:)=t(3,:)/124.7; T(4,:)=t(4,:)/89.4; pr=[0,1;0,1;0,1;0,1]; net = newff(pr,[1 1],{'logsig' 'purelin'}); net.trainparam.show=10000; %每次循环50次 net.trainParam.epochs = 500; %最大循环500次 net.trainparam.goal=0.001; %期望目标误差最小值 net = train(net,T,P); %对网络进行反复训练 Y = sim(net,T) InputWeights=net.iw{1,1}; InputWeights ——————————————————————————————————————————————— %取隐层到输出层的权值 LayerWeights=net.lw{2,1}; LayerWeights b1=net.b{1} b2=net.b{2} %画图 figure(1) plot(P,':o') hold on plot(Y,'-*') error=sum(abs(Y-P)./P) 17 110.690.495.273.65076 附录2:第三问动态规划乘法模型 #include<cstdio> #include<iostream> using namespace std; #include<cstring> #include<algorithm> #include<cmath> #define M 474 double dp[5][1000]; int use[5][1000]; ——————————————————————————————————————————————— double a[]={148.1,111.8,124.7,89.4}; double b[]={0.2090,0.7025,0.2199,0.1523}; double val(int x,int k) { double xx=x*1.0; return pow(xx/a[k-1],b[k-1]); } double gg(int x,int y,int z,int w) { return pow(x*1.0/148.1,0.2090)*pow(y*1.0/111.8,0.7025)*pow(z*1.0/124.7,0. 2199)*pow(w*1.0/89.4,0.1523); } int main() { int e,f,g,h; memset(dp,0,sizeof(dp)); for(int i=0;i<90;i++){dp[4][i]=pow(i*1.0/a[3],b[3]);use[4][i]=i;} for(int i=3;i>=1;i--) { for(int j=1;j<=M;j++) ——————————————————————————————————————————————— { for(int v=0;v<j&&v<=a[i-1]+1;v++) { if(dp[i][j]<dp[i+1][j-v]*val(v,i)) { use[i][j]=v; dp[i][j]=dp[i+1][j-v]*val(v,i); } } } } printf("%lf\n",dp[1][M]); e=use[1][M]; 18 f=use[2][M-e]; g=use[3][M-e-f]; h=use[4][M-e-f-g]; printf("%d %d %d %d",e,f,g,h); } 附录3:第三问动态规划加法模型 #include<cstdio> #include<iostream> ——————————————————————————————————————————————— using namespace std; #include<cstring> #include<algorithm> #include<cmath> #define M 474 double dp[5][1000]; int use[5][1000]; double a[]={148.1,111.8,124.7,89.4}; double b[]={0.0525,0.5575,0.2909,0.0607}; double val(int x,int k) { double xx=x*1.0; return xx/a[k-1]*b[k-1]; } double gg(int x,int y,int z,int w) { return (x*1.0/148.1*0.0525)+(y*1.0/111.8*0.5575)+(z*1.0/124.7*0.2909)+(w* 1.0/89.4*0.0607); } int main() ——————————————————————————————————————————————— { int e,f,g,h; memset(dp,0,sizeof(dp)); for(int i=0;i<90;i++){dp[4][i]=(i*1.0/a[3]*b[3]);use[4][i]=i;} for(int i=3;i>=1;i--) { for(int j=0;j<=M;j++) { for(int v=0;v<j&&v<=a[i-1]+1;v++) { if(dp[i][j]<dp[i+1][j-v]+val(v,i)) { use[i][j]=v; dp[i][j]=dp[i+1][j-v]+val(v,i); 19 } } } } printf("%lf\n",dp[1][M]); e=use[1][M]; f=use[2][M-e]; ——————————————————————————————————————————————— g=use[3][M-e-f]; h=use[4][M-e-f-g]; printf("%d %d %d %d",e,f,g,h); } 附录4:问题四动态规划模型求解 #include<cstdio> #include<iostream> using namespace std; #include<cstring> #include<algorithm> #include<cmath> #define N 500 int M; double H=50; double ETm[10]={0,148.1,111.8,124.7,89.4}; double xishu[10]={0,0.2090,0.7025,0.2199,0.1523}; double h[10]={0,70},ETa[10]; double P[10],k,Hmax=100; int m[10],mhat[10]; double gg(double x,double y,double z,double w) { return ——————————————————————————————————————————————— pow(x*1.0/148.1,0.2090)*pow(y*1.0/111.8,0.7025)*pow(z*1.0/124.7,0. 2199)*pow(w*1.0/89.4,0.1523); } int main() { FILE * f1=fopen("1.txt","r"); for(int i=1;i<=21;i++) { fscanf(f1,"%lf%lf%lf%lf%lf",&P[1],&P[2],&P[3],&P[4], &k); printf("%d年 :\n四个时期降水量:%.0lf %.0lf %.0lf %.0lf 水利用系数:%.2lf\n",2011-i,P[1],P[2],P[3],P[4],k); double maxx=-10000; memset(ETa,0,sizeof(ETa)); memset(h,0,sizeof(h));h[1]=70; 20 memset(P,0,sizeof(P)); memset(m,0,sizeof(m)); memset(mhat,0,sizeof(mhat)); //h[10]={0,70}; ——————————————————————————————————————————————— for(int v=101;v<=474;v++) { //M=v*0.1*474; for(int m1=0;m1<=v;m1++) { for(int m2=0;m2<=v-m1;m2++) { for(int m3=0;m3<=v-m1-m2;m3++) { m[1]=m1;m[2]=m2;m[3]=m3;m[4]=v-m1-m2-m3; int flag=0; for(int i=1;i<=4;i++) { if(h[i]>Hmax)h[i]=Hmax; if(h[i]+P[i]>Hmax)P[i]=Hmax-h[i]; ETa[i]=(h[i]+m[i]+P[i])*k; h[i+1]=h[i]+P[i]+m[i]-ETa[i]; if(ETa[i]>ETm[i])flag=1; } if(flag==1)continue; //printf("%lf\n",gg(ETa[1],ETa[2],ETa[3],ETa[4])-M*1.0/N); if(gg(ETa[1],ETa[2],ETa[3],ETa[4])*0.65-v*1.0/474*0.35>maxx) { ——————————————————————————————————————————————— maxx=gg(ETa[1],ETa[2],ETa[3],ETa[4])*0.65-v*1.0/474*0.35; for(int i=1;i<=4;i++)mhat[i]=m[i]; } } } } } printf("%lf\n",maxx); printf("四个时期灌水 量:%d %d %d %d\n",mhat[1],mhat[2],mhat[3],mhat[4]); printf("总灌水量:%d\n",mhat[1]+mhat[2]+mhat[3]+mhat[4]); //-M*1.0/N*0.1 } 21 } 22 ———————————————————————————————————————————————