MATHEMATICA在数学建模中的应用

2005 年4 月第2 期陈玲菊 范建生 2 M athematics 在数学建模中的应用 M athematica在数学建模中的应用 The Appl i cation on Mathematica in Mathematical Modeling 陈玲菊 1 阂江学院数学系 福建福州350108 2 范建生 2 福建商业高 等专科学 校 福建福州350012 内 容提要本文是续 数学软件M athematics及其概率统计应用 闽 江学院学报 第24 卷第五期 之后 从数学建模的角度出 发 结合实例介绍了M athem atics 在数学 建模求解过程中的应用 从而进一步阐明 数学 软件M athem atics 的多 方面的 运用 关 键 词M athematica 数学建模 费马原理 中图分类号 0245文献标识码 A文章编号 1008一 4940 2005 02 0053一 03 一 简介 M athematics是一种集数学计算 处理与分析为一身的软件 它具有的强大的数值计算 符号计算 数学图形的绘制甚至动 画制作的功能 使之在科学研究单位和教学单位中有广泛的用 途 近年来 数学建模的蓬勃发展 由于建模过程往往要依靠 计算机在计算方面的强大功能 更有效地实现大量的数据 图 像处理 而利用M athematics 编写的程序实现这些过程本身并不 复杂 而且程序执行的效率比较高 本文将就这一方面结合实 例 谈谈M athematics在数学建模中的应用 二 数学建模的含义 什么是数学模型 这方面的定义很多 概括地说 就是用 数学语言和方法对实际问题的抽象和描述 牛顿的万有引力 定律便是数学模型的一个很好的例子 学习建模就是要学会 怎样用自 己学到的数学和计算机去解决实际问题 一般来说 一个完整的数学建模过程主要由三部分组成 1 用适当的数学方法对实际问题进行描述 2 采用各种数学和计算机手段建立并求解模型 3 从实际的角度分析模型的结果 考察是否具有实际意 义 三 几个应用M athematics 求解数学模型的例子 例一 固态酒精发酵动力模型 下面利用M athematics 的数据分析函数对固态酒精发酵过 程所形成产物的规律进行研究 为简单起见只举确定酒精浓度 与发酵时间的关系这一部分 以 温度为不变参数 取 30 接种量 10 糖化酶用量 2801 U 克为原料 实验纪录如下 表 1 t 的关系为 a bt ct2 1 时 之 I L X 利用M athematics编程计算得到未知数 b c a 日 的估计 值 从而 16一 11t 11产 1 0 00408342 0 163047e 2 r 尹 1 L l 一 X 方程 2 的图形见codel 中的Out 4 根据所建立的方程 2 我们可确定何时酒精浓度取最大 为实现发酵过程的自 动 控制提供理论支持 codel 城l 二values 二 1251 1 6 134E 1 7 1381 1 istplot values II 1 16 8 140 2 38 9 141 13 821 14 112 15 1 4 0 1 2 0 1 0 0 8 0 印 4 0 Out l In 2 4 口Graphics口 二 po ints二 1 1 61 1 2 381 1 3 821 时间123456789 测量值 163882112125134138140141 可作出如下分析 1 由 表1画出散点图 根据codel 中的out I 我们可假设酒精浓度x与发酵时间 fl二 FindFit points a bt cY 2 Ia b cI t Out 2 二 a 16 b 一 11 c ll In 3 二 points 11 4 11 21 1 5 1 251 6 1341 1 7 1 381 1 8 1 401 1 9 141 I 二 FindFit points l a RF xp 一 t Ia R 1 t out 3 引a 0 00408324 0 0 1630471 收稿日 期二 2004 12 3 作者简介 陈玲菊 1976 女 闽江学院 学系助教 范建生 191 7 4 男福建商业高等专科学校讲师 54福建商业高等专科学校学报2005 年4 月 1 n 4 PI Plot 16一 llt 11Y 2 It 1 31 p2二 Plot 1 0 00408324 0 163047Exp 一 t 川 t 3 91 Show pl p2 2 0 0 1 5 0 1 0 0 5 0 问题一要求在满足一定设计规范下最节省功率的线光源 长度 我们认为这是要求在满足B点强度不小于额定值的二 倍和C点强度不小于额定值的条件下求功率的最小值 即将 问题一转化为有约束条件下的求极值的数学问题 3 模型建立卜 Out 4 口Graphics口 在这个模型中 利用M athematica 中的Listplot 画出由散点图 在坐标系中的图像 然后再用FindFit分别对分段函数进行参数 估计 分别求出a b 及二 日 的 值 可以 看出 这程序过程比 较 简单 很容易实现 下面一个例子取自2003年的全国 建模比 赛中一题 由于篇 幅所限 只对其中一部分进行讨论 例二 车灯光线源的优化设计汽车前灯为一旋转抛物面 线光源经过车灯的焦点F 在对称轴垂直的水平方向对称放置 分布均匀 测试车灯反射光的设计规范为 1 测试环境 测试屏置于F 点正前方25m 的A点 垂直于 FA 在屏上过A作平行于地面的直线 并在直线上A的同侧取 点B C 使A C 2AB二 2 6m 2 达标条件 只考虑一次反射 当C点的光强度不小于 某一额定值T 则 B点的光强度不小于2T 我们的主要工作是 1 在设计规范的约束下 给出线光源功率最小时的长度Lo 2 在有标尺的坐标系中 就1 的结果 在测试屏上画出反 射光的亮区 3 讨论给定的设计规范 并给出我们的改进方案和建议 1 基本假设 1 假设车灯的光是单色光 并且本题目 不考虑干涉衍射的 情况 2 物面镜是理想的镜面反射 不考虑漫反射情况且不吸收 能量 线光源是理想的线光源 并且发光强度均匀分布 3 只考虑一次反射 忽略二次反射和多次反射情况 2 题目 分析 车灯在测试屏上产生的照明可以分为两部分 由反射产生 的照明和由直射光产生的照明 为了计算由线光源的照明强 度 首先必须找到光源发出点 反射抛物面以及光屏上被照射 点的几何关系 我们分别从费马原理出发 建立模型描述这种 几何关系 称为光程极值模型 为了描述照射屏上的能量情况 我们必须考察光线的能流 密度 平行光和球面光的能流密度表达形式形式简单 且为我 们熟知 将线光源分为若干小块 认为每小块近似成为点光源 考虑 点光源发出的光为球面光 显然直射光照明产生的能量 可以 用球面波求解 点光源发出的波经抛物面反射后不再是 球面光 但是经平面反射的光仍然是球面光 我们可以把抛物 面分成很小的面元 这时每个面元近似认为是平面 可以按照 球面波计算 实际中这种处理法非常简洁 而且尤其在远轴 情况下非常符合精确解 图 1 如图1 以 抛物面的顶点为原点 与线光源平行的方向为x 方向 抛物面的轴的方向为 方向建立坐标系 设抛物面的方 程为 X 尹 cz 其中 的 y 为 抛 物 面 上 点 的 坐 标 值 为 常 数 并 且 有 c 二4f f为抛物面的焦距 另外假设 光屏到抛物面顶点的距 离为Z 在 线 光 源 上 坐 标为p m 0 f 发出 的 光 经 过 抛 物 面 上P 2 x y z 点 反 射 后 照 射到 光 屏上坐 标为P n 0 1 的 点 由 于 P 的y坐 标 是0 即P a将 被局限 在 一 平面 与 光 屏平面的 交 线 上 所以我们所建立的模型是一维的 4 求解的目 标 对 于 给 定的 光源 点p m 0 f 和 光 屏 上的 接受 点P 2 x y z 求出 反 射 点P a的 坐 标值 并求得单位长 度的 线光源通过P Z点 的 反 射 对P 点 处 光 强的 贡 献 然 后通 过 对 整个 线光源 进 行积 分以求得光屏上的过 中心点 的水平线上的某一点的光的相 对强度值 从而求得整条直线上的光强分布情况 但在这里 我们只利用分析方法 看看能走多远 如果读者有兴趣可以找 有关建模比 赛的材料 初步的解析求解 我们的出发点是几何光学中的费马原理 费马原理是这样 描述的 参考书目 1 的第20页 实际存在的光的路径应当是 使得这条路径取得极值的路径 极大值 极小值 均值都有可 能 光程可以 表示为 len V 一 x 2 尹 f Z 2 V 一 Y 一 1 2 在数学上 利用拉各朗日 乘数法 费马原理等效于 at ax二 at a y二 at az二 其中 t二 len k x 尹 一 cz 括号里的式子是由于如下的约束条件 2005年4 月第2 期陈玲菊 范建生 2 M athematica 在数学建模中的应用 3 尹二 cz 满足三 个 偏 微分 方程的 点 x y z 即 为 实 际 的 光反 射 点 由 a t 8 x二 得 到 y 2 4 不雨1 m x 2 y f z 2 去除虚数根 将越界的根除去 因为 由抛物面的尺寸必须有 一 仍 x 03 通过这三个步骤 我们已经得到了实际的反光点的坐标 以下程序就是要解上述方程组 求x的过程 签釜关爷釜补餐关爷釜共爷苦爷赞 爷苦书爷资 一今一昌乌一一一 0 V n一 x 2 尹 l 一 Z 可以 证明在题目 所给定的参数下这个方程得到 y 0 将这个式子代人拉各朗日 乘数法得到的方程组可知 n 一 x m一 2十 f一 爷 2 n一 x 2 f 3 n一 x 一 二 二 二 二 二 二 二 二 二 二 二 二 二 二 2x f 一二 c X 1一二 一 一 一 c 二 0 2 m一 x 2 f一 x 2 C n一 x 2子 l 2 C 这里我们用到了R 二 这个结论 将上式移项 去分母可 以得到 n一 1 一 2 x 2 一 m 2 c 一 2 c fx 2 2 V 一 z f一R 2c 一cn c一2 c lx 2 2 到此为止 我们按照费马原理作等效推导 没有出现增根 或者漏根 我们对上面的方程移项平方 并且整理可得一个二次多项 式乘以一个五次多项式的积等于0 那个二次多项式是可以用 求根公式求得其根 而那个五次多项式则无法进一步分解 到 此为止 解析的方法已经结束 我们已经得到了一个相当出色 的中间结论 我们只要求解那个五次多项式即可 尽管在刚才 移项平方的 那一步会出 现增根 5 进一步的数值解 对这个模型的进一步处理就是求解那个多项式的数值解 为此 我们需要把m n l c f的数值带到那个五次多项式中以 便使系数成为定的数值 一般来说 五次多项式会有5个根 加 上那个二次多项式的两个根 我们得到7个根 也就是说 实际 的反射点的坐标值是在这7个数中 我们需要对这7个值进行 过滤 以便得到实际的反射点 我们分析出现增根的原因 进行了三步的 过滤 去除两边平方产生的增根 产生增根的方程 一 m 一2 c fx 对 A n一x 2十 1一 兰 2二 一 n Y一c C二 一2 十 A m一 2 f一 兰 2二 我们并没有将解得的x值直接代人上式看是否为零 而是 判断两边根式之前的系数的符号 这样可以有效的避免误差带 来的对判断增根的干扰 code2 拉格朗日 乘数法推导的程序 C lear m n l x y f c Len二 m一 x 2 y2 f一 z 2 V n一 x 0 2 件 z一 1 2 t ten lamda x0 2 归一 c z rl 二 D t x r2二 D t y r3二 D t z 衅二 rl y 0 r5 r3 y o z 6 Solve x 50 0 lam daf 0 二z4 z6 沼二刃 z 2 c f l 二 r8 1 r10 Factor f l r1l 二 10 r10 3 r10 2 r10 1 r1 0 2 r1 0 3 r1 2 Coef f i cient r11 1 r10 3111 r1 3二 Coef f i cient r11 1 r10 2111 r1 4二 Expand rl2 2 r10 311 21 r15二 Expand rl3 2 r10 2 2 r16 r14一 r1 5 c 4f r17二 Factor r16 r18二 Solve r170