数学模型与数学建模实验一.docx

实验报告一 学院名称：理学院 专业年级： 姓 名： 学 号： 课 程：数学模型与数学建模 报告日期：2015 年 11 月 10 日 一、实验题目 例 1.1.1 人体的体重与身高利用下表中的一组实际测量的数据，分析人体体重 W 与身高 L 的关系。 W/kg 12 17 22 35 48 54 66 75 L/m 0.86 1.08 1.16 1.35 1.55 1.67 1.78 1.85 例：1.1.4 圆珠体积利用累次积分、示性函数积分、蒙特卡洛模拟方法计算一个直径为 2cm 的球状，带有 穿过球心，直径为 0.2cm 的圆柱体的圆珠的体积，并在达到相同精度计算结果下比较这几 种方法的计算时间。 二、实验目的 1、了解 MATLAB 的基本功能，配合介绍数学建模的最基本的要素，进行误差分析，选择合 适的数学模型。 2、为了更有效的运用 MATLAB 指令，简单运用 MATLAB 编程方法，以四种不同的方式计算 圆珠体积，通过比较不同程序的计算用时，通过比较，得到计算用时最少，精度能得到要 求的算法。 三、问题陈述 1、人体的身高与体重利用下表中的一组实际测量的数据，分析人体体重 W 与身高 L 的关系。 W/kg 12 17 22 35 48 54 66 75 L/m 0.86 1.08 1.16 1.35 1.55 1.67 1.78 1.85 2、圆珠体积的计算利用累次积分、示性函数积分、蒙特卡洛模拟方法计算一个直径为 2cm 的球状，带有 穿过球心，直径为 0.2cm 的圆柱体的圆珠的体积，并在达到相同精度计算结果下比较这几 种方法的计算时间。 四、模型及求解结果 1、人体的体重与身高解：根据质量与体积、体积与长度之间的量纲关系，可以建立轮廓模型3 kL W  关键在于参数 的估计。利用实际测量的数据估计参数 ，有两种方法。 k k 估计方法 1：用各观测数据的比值 的平均数 来估计 3 i i i L W k    i i k n k 1 k 估计方法 2：用观测数据的均值的比值估计 ： k3 1 1          i i i i L n W n k 两种参数估值 不同，估计的好坏可以通过误差平方和的大小进行比较，计算结果通过图 k 像直观地观察，结果如图：（a）估计方法 1从数值结果比较，第一种参数估计方法的误差平方和为 Q1=342.5765，小于第二种参数估 计方法的误差平方和 Q2=814.3192，因此结果更好一些，但是它仍然不是使得误差平方和 最小的估计，运用拟合模型可以重新 这个问题。 2、圆珠体积 将求体积的三重积分转为累次积分可以直接得到圆珠的体积            124 . 4 2 1 1 . 0 99 . 0 99 . 0 2 0 rdr dz d dv V z       . 01 . 0 ; 1 ; , , 2 2 2 2 2        y x z y x z y x 解法 1：换成极坐标积分，采用符号函数积分，运行程序在下面的代码中，计算用了 61.083812 seconds 解法 2：用累次积分公式，采用数值函数定积分指令。因为定积分只能计算积分上下限为 常数的定积分，虽然可以计算含参量的一重积分，但是再累次积分就只能利用数值梯形积 分公式近似计算。通过代码运行可得：计算用时 116.902487 seconds 解法 3：采用示性函数积分方法，集合 的示性函数 定义为：  n x           x x x x , 0 , 1 因此，很容易采用矩阵的逻辑运算方式定义区域 的实行函数 ，又因为该圆珠可以放   x 在立方体 内，由 在立方体 上的积分就可以得到 的体积。计算用时   3 1 , 1   x   3 1 , 1   65.264352 seconds.（b）估计方法 2解法 4：利用蒙特卡洛积分，在区间 内按均匀分布规则书 n 个点，由落在区域   3 1 , 1  内点的个数 s，可以近似计算 在 中所占的比例 ，由 得体积   3 1 , 1       3 1 , 1  n s   3 1 , 1  8，于是可以得到 的近似值 ，计算用时为 31.230011 seconds.  n s 8 五、程序代码 1、人的升高和体重 >> W=[12,17,22,35,48,54,66,75]; >> L=[0.86,1.08,1.16,1.35,1.55,1.67,1.78,1.85]; >> k1=mean(W./L.^3)%估计方法 1. k1 =13.5892 >> k2=mean(W)/(mean(L))^3%估计方法 2. k2 =14.5929 >> plot(L,W, o ,L,k1*L.^3, -k , LineWidth ,2); >> figure(2) >> xlabel( L ) >> ylabel( W ) >> plot(L,W, o ,L,k2*L.^3, -k , Linewidth ,2); >> Q1=(W-k1*L.^3)*(W-k1*L.^3) %误差平方和 Q1=342.5765 Q1 =342.5765 >> Q2=(W-k2*L.^3)*(W-k2*L.^3) %误差平方和 Q1=814.3192. Q2 =814.3192 2、圆珠体积 解法 1：换成极坐标积分，采用符号函数积分 >> tic >> int(int(sym( 2*r*sqrt(1-r^2) ),0.1,1),0,2*3.14) ans =(15543*11^(1/2))/12500 >> vpa(ans,6) ans =4.12402 >> toc Elapsed time is 61.083812 seconds. 解法 2：用累次积分公式 >> tic >> f=@(z)quad(@(r)r,0.1,sqrt(1-z.^2)); >> n=40; >> x=linspace(0,sqrt(0.99),n);>> for j=1:nf1(j)=4*3.14*f(x(j));end >> v1=trapz(x,f1)%梯形积分公式 v1 =4.1233 >> toc Elapsed time is 116.902487 seconds. 解法 3：采用示性函数积分方法 >> tic >> f2=@(x,y,z)((x.^2+y.^2+z.^2)=0.01);%示性函数 >> v2=triplequad(f2,-1,1,-1,1,-1,1) v2 =4.1231 >> toc Elapsed time is 65.264352 seconds. 解法 4：采用蒙特卡洛积分方法： >> tic >>A=unifrnd(-1,1,3,10000000); >