1、1实验一 空间曲线与曲面的绘制本实验的目的是利用数学软件 Mathematica 绘制三维图形来观察空间曲线和空间曲面图形的特点,以加强几何的直观性。1、 空间曲线的绘制绘制空间曲线时一般使用曲线的参数方程,利用命令“ParametricPlot3D” 。如画出参数方程 所确定的空间曲线的命令格式为:21 ,)()(ttzyxParametricPlot3Dxt,yt,zt,t,tmin,tmax,选项例 1 画出旋转抛物面 与上半球面 交线的图形。2yx21yxz解:它们的交线为平面 上的圆 ,化为参数方程为 ,1z22 ,0 ,1sincotz下面的 mathematica 命令就是作出它
2、们的交线并把它存在变量 p 中:p、 ParametricPlot3D、Cos、t、,Sin、t、,1、,、t,0,2、 Pi、运行即得曲线如图 1 所示。在这里说明一点,要作空间曲线的图形,必须先求出该曲线的参数方程。如果曲线为一般式 ,其在 面上的0),(zyxGFxOy投影柱面的准线方程为 ,可先将 化为参),H),(H数方程 ,再代入 或 解出 )(tyx0,(zyx0,zyxF即可。 tz2、 空间曲面的绘制作一般式方程 所确定的曲面图形的 Mathematica 命令为:),(yxfzPlot3Dfx,y,x,xmin,xmax,y,ymin,ymax,选项 作参数方程 所确定的曲
3、面图形的,)( maxinmaxinvuvzyMathematica 命令为:ParametricPlot3Dxu,v,yu,v,zu,v,u,umin,umax,v,vmin,vmax,选项例 2 作出上半球面 的图形。21yxz -1-0.5 0 0.5 1-1-0.500.5100.511.52 、2解:首先我们选取绘图区间 作图,输入下面语句:1,1yxPlot3D、1、1、 x2、 y2,、x,、1,1、,、y,、1,1、运行后得到了该曲面的图形(图 6-2) ,但是在图形的前面出现了一些蓝色字体报错信息,而且图形不完整,这是因为函数 在范围 内的一些点处无定义。z,x为避免上述问题
4、,可用下面两种方法:(1) 定义一个分区域函数 ,将无定义的点赋予函数值 1:)(xf,,1,),( 22yxyxf作出该函数的图形只要键入命令:f、x_,y_、:、 If、x2、 y2、 1,1、1、 x2、 y2,1、Plot3D、f、x,y、,、x,、1,1、,、y,、1,1、运行后得图 3,可以看到该图形比上半球面多了一部分曲面的图形(即 平面上的z部分) 。但是图形比较粗糙,我们可以提高采样点数,例如取采样点数为 30,即运行命令Plot3D、f、x,y、,、x,、1,1、,、y,、1,1、,PlotPoints、 30、可得图形 4,由此可见图形已经比较精细了。-1 -0.50 0
5、.51-1-0.500.5111.251.51.752-1 -0.50 0.51-1-0.500.5111.251.51.752-1 -0.50 0.51-1-0.500.5111.251.51.752图 2 图 3 图 4(2) 采用参数方程,选取参数的范围使得区域内的每一点都有定义。对于题目中的球面有参数方程 ,我们输入命令:2,0 , ,cos1invuvzyxParametricPlot3D、Cos、u、 Sin、v、,Sin、u、 Sin、v、,1、 Cos、v、,、u,0,2、 Pi、,、v,0,Pi2、,PlotPoints、 30、运行后得图形 5。我们还可以改变参数的范围画出
6、上半球面的 部分(如图 6):4ParametricPlot3D、Cos、u、 Sin、v、,Sin、u、 Sin、v、,1、 Cos、v、,、u,0,3、 Pi、2、,、v,0,Pi、2、,PlotPoints、 30、3-1 -0.50 0.51-1-0.500.5111.251.51.752-1-0.50 0.51-1-0.500.5111.251.51.752图 5 图 63、 空间图形的叠加与平面图形类似,空间的立体图形同样可用“Show”命令,把不同的图形(曲线或曲面)叠加并在一个坐标系中显示出来。例 3 画出由旋转抛物面 与上半球面 相交所围成的立体几2yxz21yxz何图形。解
7、:这是一个组合图形。一般地,直接画出两者的图形再组合在一起。但是,这里所要的图形仅仅是两个曲面图形的一部分,因此需要有选择地画出两曲面的相应部分再组合。由于它们的交线为 ,故相应的曲面部分的参数方程为:12zyx与 。,02 , ,sinco2rtrzyx 2,0 , ,cos1invuvzyx输入以下 Mathematica 语句:t1、 ParametricPlot3D、r、 Cos、t、,r、 Sin、t、,r2、,、t,0,2、 Pi、,、r,0,1、,PlotPoints、 30、;t2、 ParametricPlot3D、Cos、u、 Sin、v、,Sin、u、 Sin、v、,1、
8、 Cos、v、,、u,0,2、 Pi、,、v,0,Pi2、,PlotPoints、 30、;Show、t1,t2、运行后即得旋转抛物面、上半球面及叠加曲面的图形(图 7) 。-1 -0.50 0.51-1-0.500.5100.250.50.751-1 -0.50 0.51-1-0.500.5111.251.51.752-1-0.50 0.51-1-0.500.5100.511.52图 7例 4 绘制由曲面 与 所围成的立体区域。22xyxz、z解:输入命令:4s1、 ParametricPlot3D、u,v,u2、 v2、,、u,、1,1、,、v,、1,2、,PlotRange、0,2、,A
9、xesLabel、“X“,“Y“,“Z“、,DisplayFunction、 Identity、;s2、 ParametricPlot3D、u2,u,v、,、u,、1,1、,、v,0,2、,AxesLabel、“X“,“Y“,“Z“、,DisplayFunction、 Identity、;s3、 ParametricPlot3D、1,u,v、,、u,、1,1、,、v,0,2、,AxesLabel、“X“,“Y“,“Z“、,DisplayFunction、 Identity、;s4、 ParametricPlot3D、u,v,0、,、u,、1,1、,、v,、1,1、,AxesLabel、“X“,
10、“Y“,“Z“、,DisplayFunction、 Identity、;Show、s1,s2,s3,s4,DisplayFunction、 $DisplayFunction、 在上述语句中,选项“DisplayFunction Identity”表示不显示图形,而“DisplayFunction$DisplayFunction”则表示显示图形。运行结果如图 8。4用动画来演示产生旋转曲面的过程。例 5 用动画演示由曲线 绕 轴旋转产生旋转曲面的过程。,0sinzyz解:该曲线绕绕 轴旋转产生的曲面方程为 ,其参数方程为z zyx22sin,输入以下命令,就可得到连续变化的 20 幅图形:,0,
11、 ,sincouzyxm、 20;For、i、 1,i、 m,i、,ParametricPlot3D、Sin、z、 Cos、u、,Sin、z、 Sin、u、,z、,、z,0,Pi、,、u,0,2、Pi、 i、m、,AspectRatio、 1,AxesLabel、“X“,“Y“,“Z“、,PlotPoints、 30、运行后得到 20 幅曲面的图形,图 8 中列举了其中的三幅。大家还可以进行动画演示,观察到旋转曲面产生的过程。0 0.25 0.5 0.75 1X00.20.40.60.8Y0123Z-1 -0.5 0 0.51X00.250.50.751Y0123Z-1 -0.5 0 0.51
12、X-1-0.500.51Y0123Z图8 -1-0.500.51X-1012Y00.511.52Z 8 、5实验习题1、 作出各种标准二次曲面的图形。2、 利用参数方程作图,作出由下列曲面所围成的立体:(1) 及 面xyxz22、Oy(2) 及01y、z3、观察二次曲面族 的图形。特别注意确定 的这样一些值,当 经过kxyz2 kk这些值时,曲面从一种类型变成了另一种类型。实验二 无穷级数与函数逼近本实验的目的是用 Mathematica 显示级数部分和的变化趋势;学会如何利用幂级数的部分和对函数的逼近以及进行函数值的近似计算;展示傅里叶级数对周期函数的逼近情况。1、 级数部分和的变化趋势在实
13、验一中,我们通过图形能够清楚地显示极限的变化趋势,而级数的和就是部分和序列的极限,下面我们采用散点图来观察级数部分和序列的变化趋势。例 1 观察级数 的部分和序列的变化趋势,并求和。12)(n解:(1)可以用以下两种方法,从图形来观察级数的敛散性。(a)利用“Table”命令生成部分和数列的数据点集后作点图,输入语句如下:s、n_、:、 Sum、1、k、 1、2、k、 1 ,、k,1,n、;data、 Table、s、n、,、n,1,400、;ListPlot、data、100 200 300 4000.7780.7820.7840.7860.7880.790.7920.784 0.785 0
14、.786 0.787 0.788-0.20.20.40.60.811.2图 1 图 2运行后见图 1。从图中可以看到级数收敛,级数和大约为 0.786。(b)将级数的所有部分和用竖直线段画出,得到类似条形码的图形,通过这种图形来看出6级数的收敛性。输入命令如下:sn、 0;n、 1;h、;m、 3;While、1n、 10、m、,sn、 sn、 、1、n、 1、2、n、 1 ;h、Append、h,Graphics、RGBColor、Abs、Sin、n、,0,1、n、,Line、sn,0、,、sn,1、;n、 ;Show、h,PlotRange、0.2,1.3、,Axes、 True、 运行后
15、见图 2。从图中可以看出,级数的和在 0.785 与 0.786 之间,并且可以通过改变m 的值来提高观察到的和的精度。(2)求和。对于这个级数,可以通过基本输入模板利用求和符号来直接求和,只要输入命令、n、1、 、1、n、 1、2、n、 1 ,运行后得到级数和的精确值: 。但并不是所有的级数都能如此求4得和,但可以利用下面的命令来求和的近似值:、n、1、 、1、n、 1、2、n、 1 、NNSum、1、n、 1、2、n、 1 ,、n,Infinity、运行后得结果均为:0.785398。为了得到更高的精度,还可以用“N”命令来求函数“Sum”所得的值。大家可以比较一下下面两条命令输出结果的差
16、异:N、Sum、1、n、 1、2、n、 1 ,、n,Infinity、,30、N、NSum、1、n、 1、2、n、 1 ,、n,Infinity、,30、2、 函数的幂级数展开如果函数 在 邻域内具有任意阶导数,则函数可以展开为 处的幂级数)(xf0 0x, 100)(10 )!)(nnn xfxaf称之为泰勒级数。特别地,当 时,称为麦克劳林级数。0例 2 将函数 展开为 的幂级数,并利用图形考察幂级数的部分和逼近函mxf)(数的情况。解:根据幂级数的展开公式,若 能展开成 的幂级数,其展示为)(fx,因此首先定义函数,再计算 点的 阶导数,最后构成和式。1)(!0)(nnxfxf 0n不妨
17、设 ,输入命令如下:m7m、 、1;f、x_、:、1、 x、m;x0、 0;g、n_,x0_、:、 D、f、x、,、x,n、.x、 x0;s、n_,x_、:、 Sum、g、k,x0、k、 、 x、 x0、k,、k,0,n、;命令中的函数 sn_,x_表示的是函数 在 处的 n 阶泰勒级数。下面我们通过以下)(f命令观察幂级数的部分和逼近函数的情况:t、 Table、s、n,x、,、n,20、;p1、 Plot、Evaluate、t、,、x,、1、2,1、2、;p2、 Plot、1、 x、m,、x,、1、2,1、2、,PlotStyle、 RGBColor、0,0,1、;Show、p1,p2、运
18、行结果如图 3。从图中可以看到,当 n 越大时,幂级数越逼近函数。-0.4 -0.2 0.2 0.40.60.81.21.41.61.82-0.4 -0.2 0.2 0.40.81.21.41.61.82-0.4 -0.2 0.2 0.40.60.81.21.41.61.82图 33、傅里叶级数设 是以 2T 为周期的周期函数,在任一周期内, 除在有限个第一类间断点)(xf )(xf外都连续,并且只有有限个极值点,则 可以展开为傅里叶级数:)(xf,其中 ,10 )sincos(2n TbxaTn ndxfba,321 ,si)(10co且傅里叶级数在任一点 处收敛于 。0x20)(0xxf下
19、面通过例子来研究傅里叶级数的生成以及级数部分和逼近函数的情况。例 3 画图观察周期为 、振幅为 1 的方波函数 展成的傅里2xxf0,1)(叶级数的部分和逼近 的情况。)(xf解:根据傅里叶系数公式可得: ,0)(10dxfa, ,cos)cos(10ndxnan sin)sin(10xdxbn故输入以下命令,从输出的图形(图 4)观察傅里叶级数的部分和逼近 的情况:(f8f、x_、:、 Which、2、Pi、 x、 、Pi,1,、Pi、 x、 0,、1,0、 x、 Pi,1,Pi、 x、 2Pi,、1、;a、n_、:、 Integrate、Cos、nx、,、x,、Pi,0、Integrate
20、、Cos、nx、,、x,0,Pi、Pi;b、n_、:、 Integrate、Sin、nx、,、x,、Pi,0、Integrate、Sin、nx、,、x,0,Pi、Pi;s、x_,n_、:、 a、0、2 、 Sum、a、k、 Cos、kx、 b、k、 Sin、kx、,、k,1,n、;g1、 Plot、f、x、,、x,、2、Pi,2、Pi、,PlotStyle、 RGBColor、0,0,1、,DisplayFunction、 Identity、;m、 18;For、i、 1,i、 m,i、 2,g2、 Plot、Evaluate、s、x,i、,、x,、2、Pi,2、Pi、,DisplayFunc
21、tion、 Identity、;Show、g1,g2,DisplayFunction、 $DisplayFunction、-6 -4 -2 2 4 6-1-0.50.51-6 -4 -2 2 4 6-1-0.50.51-6 -4 -2 2 4 6-1-0.50.51-6 -4 -2 2 4 6-1-0.50.51-6 -4 -2 2 4 6-1-0.50.51-6 -4 -2 2 4 6-1-0.50.51-6 -4 -2 2 4 6-1-0.50.51-6 -4 -2 2 4 6-1-0.50.51-6 -4 -2 2 4 6-1-0.50.51图 4从图表可以看出,n 越大逼近函数的效果越
22、好,还可以注意到傅里叶级数的逼近是整体性的。例 4 将函数 展开成周期为 6 的傅里叶级数。1)(4xf解:输入命令如下:9fourier、f_,T_,k_、:、 Module、a,b,i,t,s,g1,g2、,a、0、 Integrate、f,、x,、T,T、T;s、 a、0、2;For、i、 1,i、 k,i、,t、 i、 Pi、T;a、i_、:、 Integrate、f、 Cos、t、 x、,、x,、T,T、T;b、i_、:、 Integrate、f、 Sin、t、 x、,、x,、T,T、T;s、 s、 a、i、Cos、t、 x、 b、i、Sin、t、 x、;Print、s、;Plot、
23、Evaluate、s、,、x,、T,T、,DisplayFunction、 Identity、;f、 x4、 1;T、 3;n、 8;g、 Plot、f,、x,、T,T、,PlotStyle、 RGBColor、0,0,1、;For、j、 1,j、 n,j、 2,p、 fourier、f,T,j、;Show、g,p,DisplayFunction、 $DisplayFunction、在上面的 Mathematica 程序中,先定义了一个产生傅里叶级数的前 k 项部分和的函数,在后面的“For”循环中连续调用该函数,输出了 的傅里叶级数的前)(xf项部分和函数,并画出了 图形及同一坐标下 与其傅
24、里叶级数前),7531(nf )(xf7 项和函数的图形(如图 5) 。-3 -2 -1 1 2 3246810-3 -2 -1 1 2 351015202530图5实验习题81、 观察级数 的部分和序列的变化趋势,并求和1!n2、 改变例2中m及 的数值来求函数的幂级数及观察其幂级数的逼近函数的情况。0x3、 观察函数 展成的傅里叶级数的部分和逼近 的情况。xf,10)( )(xf实验三 最小二乘法10在科学研究和实际工作中,常常会遇到这样的问题:给定两个变量 x, y 的 m 组实验数据 ,如何从中找出这两个变量间的函数关系的近似解析表达),(),(,21myxyx式(也称为经验公式) ,
25、使得能对 x 与 y 之间的除了实验数据外的对应情况作出某种判断。这样的问题一般可以分为两类:一类是对要对 x 与 y 之间所存在的对应规律一无所知,这时要从实验数据中找出切合实际的近似解析表达式是相当困难的,俗称这类问题为黑箱问题;另一类是依据对问题所作的分析,通过数学建模或者通过整理归纳实验数据,能够判定出 x 与 y 之间满足或大体上满足某种类型的函数关系式 ,其中),(axf是 n 个待定的参数,这些参数的值可以通过 m 组实验数据来确定(一),(21aa般要求 ) ,这类问题称为灰箱问题。解决灰箱问题的原则通常是使拟合函数在 处m ix的值与实验数值的偏差平方和最小,即 取得最小值。
26、这种在方差意义ni iiyaxf12),(下对实验数据实现最佳拟合的方法称为“最小二乘法” , 称为最小二乘解,na,21称为拟合函数。下面对于“最小二乘法”通过两个例子来说明。),(axfy一、线性拟合问题例 1在某化工厂生产过程中,为研究温度 x(单位:摄氏度)对收率(产量)y (%)的影响,可测得一组数据如下表所示,试根据这些数据建立以 x 与 y 之间的拟合函数。温度 x 100 110 120 130 140 150 160 170 180 190收率 y(%) 45 51 54 61 66 70 74 78 85 89解:(1)确定函数的类型为此,应对给定的数据以 x 为横坐标,y
27、 为纵坐标作出散点图。输入如下语句:x、 Table、100、 10、 i,、i,0,9、;y、45,51,54,61,66,70,74,78,85,89、;xy、 Table、x、i、,y、i、,、i,1,10、ListPlot、xy,PlotStyle、 PointSize、0.015、运行后可得到数据表和图 1。从图 1 种可以看出这些点近似地落在一条直线周围,可以认为在 y 和 x 之间存在着线性关系,之所以不完全落在直线上,是因为观察数据本身存在的误差。下面我们就用最小二乘法求出与这些数据点最接近的直线方程。11、100,45、,、110,51、,、120,54、,、130,61、,
28、、140,66、,、150,70、,、160,74、,、170,78、,、180,85、,、190,89、 120 140 160 18060708090图 1(2)利用 Mathmatica 求最小二乘解设直线方程式为 ,其中 a, b 是待定系数,于是,可以把拟合函数在 处的axy ix值与观察值的偏差平方和表示为关于待定系数 a, b 的二元函数:,ni iiyxQ12)(),(按照二元函数求极值的理论,其最小值应满足方程组:,即 。0)1()(21ni iii iiybaxbQx0)(1ni iii iiybax于是,输入以下 Mathematica 语句便可求得待定系数 a, b:q
29、、a_,b_、:、 Sum、ax、i、 b、 y、i、2,、i,1,10、Solve、D、q、a,b、,a、 0,D、q、a,b、,b、 0、,、a,b、运行后可得解:也可以用“NSolve”命令得到 a, b 的小数解: 、a、 0.48303,b、 、2.73939、。这样,我们得到方程 ,这就是前面所说的经验公式。739.2480.xy在此例题中,观察的数据点大致呈一直线状,如果在实际问题中数据点落在某一条曲线周围,则要根据曲线的形状来确定拟合函数的类型。比如考虑用 x 的 n 次多项式来拟合,这成为多项式拟合。此时仍可用最小二乘法,nn xaxaP210)(考虑 元函数 取最小的必要条
30、件,令此函数对i iinnyxPQ12210 )(),(各个参数的偏导等于 0,解一个 元的方程组便可求得这些参数的最小二乘解。二、可化为线性拟合的曲线拟合问题在许多场合下,拟合函数不具有线性形式,但是由实际经验或相关的学科理论,能够提供拟合函数的可取类型,而且可以通过适当的变量代换将拟合函数线性化,同样可以建立经验公式。(1)模型 可以用变量替换 将函数化为线性函数:bxaeyxXyY,ln、a、 7971650,b、 、 452165、12。bXaYln(2)模型 可以用变量替换 将函数化为线性函数。xylnxXyYln,(3)模型 可以用变量替换 将其化为 z 和 x 之间的线性拟合。)
31、()(z例 2研究黏虫的生长过程,可测的一组数据下表所示。温度 t 11.8 14.7 15.4 16.5 17.1 18.3 19.8 20.3历期 N 30.4 15 13.8 12.7 10.7 7.5 6.8 5.7其中历期 N 是指卵块孵化成幼虫的天数。昆虫学家认为在 N 与 t 之间有关系式:,其中 k, c 为常数。试求最小二乘解。t解:作变换 ,由此把 N, t 的关系式化为了关于 x, y 的线性关kcbatxy,1系式: 。与例 1 相同,输入以下语句求解参数 a, b:xktt、11.8,14.7,15.4,16.5,17.1,18.1,19.8,20.3、;w、30.4
32、,15.0,13.8,12.7,10.7,7.5,6.8,5.7、;y、 1、w;x、 t;xy、 Table、x、i、,y、i、,、i,1,8、;q、a_,b_、:、 Sum、ax、i、 b、 y、i、2,、i,1,8、Solve、D、q、a,b、,a、 0,D、q、a,b、,b、 0、,、a,b、运行后的结果为: 、a、 0.016481,b、 、0.175433、。最后为了求出参数 k, c 以及经验公式,还需要回代变量,为此输入命令:A、a,b、.%;k、 1、A、1,1、c、 、k、 A、1,2、运行后得到参数 k, c 的值为:k =60.6758; :c=10.6445。最后,拟
33、合曲线方程为:。645.1078tN为了比较得到的拟合函数和已知的数据点,我们再在同一坐标下绘出数据点的散点图及拟合函数的图形,输入语句为:13data、 Table、t、i、,w、i、,、i,1,8、;t1、 ListPlot、data,PlotStyle、 PointSize、0.02、,DisplayFunction、 Identity、;f、x_、:、 k、x、 c、;t2、 Plot、f、x、,、x,11,25、,AxesOrigin、11,0、,DisplayFunction、 Identity、;Show、t1,t2,DisplayFunction、$DisplayFunctio
34、n、图 2运行结果如图 2 所示。实验习题1为测定刀具的磨损速度,每隔一小时测量一次刀具的厚度,由此得到以下数据:时间 t 0 1 2 3 4 5 6 7厚度 y 27.0 26.8 26.5 26.3 26.1 25.7 25.3 24.8试根据这组数据建立 y 与 t 之间的拟合函数。2一种合金在某种添加剂的不同浓度下进行实验,得到如下数据:浓度 x 10.0 15.0 20.0 25.0 30.0抗压强度 y 27.0 26.8 26.5 26.3 26.1已知函数 y 与 x 的关系适合模型: ,试用最小二乘法确定系数 a,b,c,并2cxbay求出拟合曲线。3在研究化学反应速度时,得到下列数据: ix3 6 9 12 15 18 21 24iy57.6 41.9 31.0 22.7 16.6 12.2 8.9 6.5其中 表示实验中作记录的时间, 表示在相应时刻反应混合物中物质的量,试根据这些ixiy数据建立经验公式。14 16 18 20 22 24255075100125150