1、Maple 简介,周义仓 西安交通大学理学院应用数学系Tel. 82668741(O)中2-2210Email: ,主要内容和参考文献1.1 主要内容参考资料Maple 软件包介绍 引言,数值计算 ,代数运算 , 图形 ,微积分,线 性 代 数 , 微分方程, Maple 编程 Maple 应用举例 最速降线问题,人体含铅量,1.2 参考资料Maple与微分方程(复印);周义仓、靳祯、秦军林,常微分方程及其应用,科学出版社,2003年6月。精英科技,孙非,Maple6实例教程,中国电力出版社,2001。李强,Maple8基础教程,中国水利水电出版社,2004David Barrow等, Sol
2、ving Differential Equations with Maple V, Brooks/Cole Publishing Company, 1998.http:/202.112.84.70/resource/maple.htmhttp:/202.204.240.30/maple/tour/index.htm,Maple9教程“四部曲” Maple 9 Getting Started Guide , 834KB http:/bolide.lamost.org/sheaweb/files/m9GettingStartedGuide.pdfMaple 9 Learning Guide, 4.
3、49MB http:/bolide.lamost.org/sheaweb/files/m9LearningGuide.pdf Maple 9 Introductory Programming Guide, 7.08MB http:/bolide.lamost.org/sheaweb/files/m9IntroductoryProgrammingGuide.pdf 4. Maple 9 Advanced Programming Guide , 3.63MB http:/bolide.lamost.org/sheaweb/files/m9AdvancedProgrammingGuide.pdfMa
4、ple应用http:/ 软件包介绍2.1 引言 Maple是Waterloo大学开发的集数值计算、符号运算和图形显示于一体的软件包。它可以把人们从烦琐的代数运算中解脱出来,使得我们将主要精力集中与分析和解决问题的思路与方法中去。数学运算分为数值运算、符号运算、逻辑运算等。 数值运算: 如函数求值、方程求根、矩阵的特征值、特征向量等,数值计算是一项简单烦琐的工作,计算机的出现解决了数值计算的困难,使得大规模、高速度的计算可以实现,使得数学在天气预报、油藏模拟、航天等领域发挥了巨大的作用。,符号运算(代数运算):这是一种更加智能化的计算,所处理的符号可以代表各种数、代数式、数学结构(集合),代数运
5、算是在一些代数规则下所进行的数学处理,主要是寻求一个简单完美的公式解答。 用计算机进行符号运算是数学与计算机领域中一个新的发展方向,是人们用计算机代替大脑进行数值运算后用计算机代替大脑进行代数运算的成功实践,是数学走向机械化的重要步骤,这也使得计算机更加智能化。 目前数学软件包很多,如 MathematicaMapleMatlabMathCADSASLindo,Maple 的历史1980年Wterloo大学符号计算研究小组成立;1985年:Maple3.3商业版发行;1992年:Windows操作系统下的Maple V Release2面世;1994年,Maple V Release3发行;1
6、996年,Maple V Release4发行; 1997年,Maple V Release5发行; 1998年,Maple V Release5.2,5.2发行; 2000年,Maple 6发行; 2001年,Maple 7发行; 200?年,Maple 8发行; 200?年,Maple 9发行; 200?年,Maple 9.5发行;,Maple提供了一整套问题解决的数学环境 它支持广泛、大量的数学操作如:数值分析、符号代数、图 形图像。 Maple是一个广泛应用的符号计算软件,它拥有强大的功能,主要表现在集符号运算、数值计算、可视化和程序设计于一体,这些功能是通过Maple提供的线性代数程
7、序包、微积分方程程序包、统计程序包、偏微分方程程序包以及画图程序包等实现的。 在和别的数学软件相比,Maple在处理微分方程中问题的功能特别好。,在线帮助 Maple 允许通过如下方式获得信息 : 在由建菜单选择 Help on XXXX,或 ctrl+F1得到鼠标指针位置文 字的帮助 。 在 help菜单选择 Introduction使用帮助导航器, 获得层级结构导航,单击感兴趣的标题直到出现相应帮助。 在help菜单选择Topic Search 进行主题搜索。在工作表中可使用?命令获得明确主题的帮助页。例如 ?plot 打开plot命令 的帮助页。 在 help菜单选择Full Text
8、Search可以进行全文关键字检索。 在 help菜单选择History可以方便的回访已访问 的帮助页。,2.2 数值计算 Maple可视为功能强大的计算器。请输入下面的指令:32*1213;200!;ifactor(%);expand(%); (230)/(320)*sqrt(3);evalf(%);evalf(Pi,500);Sum(1+i)/(1+i4),i=1.100);sum(1+i)/(1+i4),i=1.100);evalf(%);Sum( 1/k2, k=1.infinity ); value(%); Product( (i2+3*i-11)/(i+3), i=0.10 );
9、value(%); evalf( % , 50 );,2.3 代数运算 展开、分解、化简表达式 expr := (x+y)15; expand(expr); factor(%); simplify(cos(x)5+sin(x)4+2*cos(x)2-2*sin(x)2-cos(2*x) ); normal( (x3-y3)/(x2+x-y-y2) ); 定义变量 expr1 := (41*x2+x+1)2*(2*x-1); expr2 := expand(expr1); eval(expr2 , x=1 ); top := expr2; bottom := expand(3*x+5)*(2*x
10、-1); answer := normal( top/bottom );,表达式变形 my_expr := (a*x2+b)/(x*(-3*x2-x+4); my_expr := (a*x2+b)/(x*(-3*x2-x+4); 定义函数 f := x - x2+1/2 ; f(2); f(a+b); 使用 unapply命令将表达式转化为函数g := unapply(x2 + 1/2, x); 解方程(组) eqn:=x3-1/2*a*x2+13/3*x2=13/6*a*x+10/3*x-5/3*a; solve( eqn, x ); 为验根我们计算方程在特殊点x的值 eval( eqn ,
11、 x=1/2*a );,定义如下方程组;eqn1 := a+2*b+3*c+4*d+5*e=41; eqn2 := 5*a+5*b+4*c+3*d+2*e=20; eqn3 := 3*b+4*c-8*d+2*e=125; eqn4 := a+b+c+d+e=9; 用变量e来表示其他未知数 a, b, c, d得到一组解 solve( eqn1, eqn2, eqn3, eqn4, a, b, c, d ); 解不等式ineq := x+y+4/(x+y) 0, x-y x*sin(a*x) + b*x2; #定义函数diff( f(x), x ); #求导 int(f(x), x); #积分i
12、nt(f(x), x=1.2); #积分expr := (2*x+3)/(7*x+5); #定义表达式Limit( expr, x=infinity ); #显示求表达式的极限limit( expr, x=infinity ); #给出表达式的极限expr := sin(4*x)*cos(x): #定义表达式approx1 := series( expr, x=0 ); #展开为级数poly1 := convert( approx1, polynom ); #转化为多项式 plot( expr, poly1, x=-1.1, y=-2.2, title = cat( convert(expr,
13、 string), vs. Series Approximation ) ); #比较展开式与原函数的图像,2.6 线 性 代 数 使用with(linalg)命令载入线 性代数工具包 restart; with(linalg); A:=matrix( 3, 3,1/2, -1/3, 2, -5, 14/3, 9, 0, 11, -5/6 );#定义矩阵det(A); #计算行列式值 inverse(A); #计算逆矩阵 B:=matrix(3,3, 1/2, 0, -2, sin(theta),1,phi2, 0,phi-1, 3/4 ); C := multiply( A, B ); #求
14、矩阵 A、B的积并存于C det(C); #计算行列式值 eigenvals(A); #特征值 eigenvects(A); #特征向量 vandermonde( s, t, u, v, w ); #Vandermonde矩阵 A := matrix( 3, 3 ); print(A);,2.7 微分方程解微分方程命令格式:dsolve(equn,y(x);dsolve(equn,conds,y(x);equn为方程conds为条件例ode1 := diff(y(x),x)-y(x)-cos(x);ans1 := dsolve(ode1,y(x);ans2 := dsolve(ode1,y(0
15、)=1,y(x);ode3 := diff(y(x),x)-y(x)2+y(x)*sin(x)-cos(x);ans3 := dsolve(ode1);,解微分方程组命令格式: dsolve(sysODE, ICs, funcs)sysODE为方程组ICs为条件组例restart;sys:= diff(x(t),t)=2*x(t)+y(t),diff(y(t),t)=3*x(t)+4*y(t);ans1:= dsolve(sys, x(t),y(t);ans2:= dsolve(sys,x(0)=0,y(0)=1, x(t),y(t);assign(%);x(t);plot(x(t),y(t)
16、,t=0.1);,第21页用逐次迭代法求下列初始值问题的解: 解. 该问题等价的积分方程为: 利用Maple去进行这些重复性的迭代:y0:=1;y1:=1+int(x2+y02,x=0.x);y2:=1+int(x2+y12,x=0.x);y3:=1+int(x2+y22,x=0.x);y4:=1+int(x2+y32,x=0.x);,第25页DEtoolsphaseportrait #画向量场及积分曲线;(diff(y(x),x)=-y(x),y(x), #定义微分方程y=-y;x=-2.2, #指出x的范围y(-2)=2,y(-2)=1,y(-2)=-2, #给出3个初始值dirgrid=
17、17,17, #定义网格点密度arrows=LINE, #定义线段类型axes=NORMAL); #定义坐标系类型DEtoolsdfieldplot(diff(y(x),x)=x2-y(x),y(x),x=-2.2,y=-2.2,dirgrid=9,9,arrows=LINE,axes=NORMAL);,第26页DEtoolsphaseportrait #画向量场及积分曲线;(diff(y(x),x)= x2-y(x),y(x), #定义微分方程x=-2.2, #指出x的范围y(-2)=1.3,y(-2)=-2, #给出3个初始值dirgrid=33,33, #定义网格点密度arrows=LI
18、NE, #定义线段类型axes=NORMAL); #定义坐标系类型,第88-91: 常微分方程的数值解原理:函数的近似展开在此基础上可以得到不同精度的算法,Euler折线及改进的Euler折线法printlev1:=0;h:=0.1;x0:=0;y0:=0.5;z0:=0.5;f1:=(x,y)-1+(y-x)2;f2:=(x,y)-2*(x-y)+2*(y-x)*(1+(y-x)2);for n from 0 to 9 dox|(n+1):=h*(n+1);y|(n+1):=y|n+h*f1(x|n,y|n);z|(n+1):=z|n+h*f1(x|n,z|n)+h2*f2(x|n,y|n)
19、/2;print (x|(n+1),y|(n+1),z|(n+1);od;,整理数据data1:=x0,y0,x1,y1,x2,y2,x3,y3,x4,y4,x5,y5,x6,y6,x7,y7,x8,y8,x9,y9,x10,y10;data2:=x0,z0,x1,z1,x2,z2,x3,z3,x4,z4,x5,z5,x6,z6,x7,z7,x8,z8,x9,z9,x10,z10;作图比较plot(data1);plot(data2);plot(data1,data2,color=blue, red);,求精确解dsolve(diff(y(x),x)=1+(y(x)-x)3,y(0)=0.5,
20、y(x);求精确解在个点的值w0:=0.5;h:=0.1;for n from 0 to 9 dox|(n+1):=h*(n+1);w|(n+1)=subs(x=x|(n+1), -1/2*(-8*x+4*x2)/(4-2*x)+1/(sqrt(4-2*x);od;整理数据data3:=x0,w0,x1,w1,x2,w2,x3,w3,x4,w4,x5,w5,x6,w6,x7,w7,x8,w8,x9,w9,x10,w10;作图比较 plot(data3); plot(data1,data2,data3,color=blue, red,green);,2.8 Maple 编程Maple提供了简单有
21、效的编程语言.2.8.1 条件语句 if if 条件 then 语句组 fiif 条件 then 语句组 else 语句组 fiif 条件 then 语句组 elif 条件then 语句组 fiif 条件 then 语句组 elif 条件 then 语句组else 语句组 fix:=10:if x2 then x2 fi; if x0 then y:=1 else y:=0 fi;,2.8.2 循环语句for,while,do for 变量名 from 初值 by 步长 to 终值 do 语句组 odwhile 条件 do 语句组 od for i from 6 by 10 to 50 do
22、print(i) od; s := 0; for i from 11 by 2 while i 0 dot:=t2+1;if t5 then break fi; od;2.8.4过程 Maple的绝大多数函数是以过程的形式存在,一个Maple语言的过程定义如下:lc := proc( s, u, t, v ) description form a linear combination of the arguments; s * u + t * vend proc;lc(3,x,4,y);,L:=proc(n:nonnegint) if n=1 then 1 elif n=2 then 1 el
23、se L(n-1)+L(n-2)fiend proc:L(1);L(2);L(3);L(4);L(5);L(6);,p22:=proc(A)p:=-(A1,1+A2,2):q:=A1,1*A2,2-A1,2*A2,1:t:=p2-4*q:s1:=is(q,positive):s2:=is(q,negative):s3:=is(t,positive):s4:=is(t,0):s5:=is(t,negative):s6:=is(p,positive):s7:=is(p,0):s8:=is(p,negative): if (s1)and(s3)and(s6)then printf(稳定结点);eli
24、f(s1)and(s3)and(s8)then printf(不稳定结点);elif(s1)and(s4)and(s6)then printf(稳定的临界或退化结点);elif(s1)and(s4)and(s8)then printf(不稳定的临界或退化结点);elif(s1)and(s5)and(s6)then printf(稳定焦点); elif(s1)and(s5)and(s8)then printf(不稳定焦点);elif(s1)and(s7)then printf(中心);elif(s2)then printf(鞍点);else printf(无法断定);end if;end;,3.
25、Maple 应用举例3.1最速降线问题 在数学史上最著名的问题之一就是最速降线问题:在给定 点 P(高 )、Q(低 )之间,找一条曲线使一质点沿该曲线无摩擦下滑 的时间最短 。建立坐标系求曲线y(x)y(x)光滑y(0)=0,y(a)=b; 下滑时间最小先计算在小段时间再积分,总时间为测试: restart;a:=1;b:=1;y:=x;g:=9.8;dy:=diff(y,x);times=int(sqrt(1+dy2)/sqrt(2*g*y),x=0.a);Time=0.63887restart;a:=1;b:=1;y:=sqrt(x);g:=9.8;dy:=diff(y,x);times=
26、int(sqrt(1+dy2)/sqrt(2*g*y),x=0.a);evalf(%);Time=0.58440,restart;a:=1;b:=1;y:=sqrt(1-x2);g:=9.8;dy:=diff(y,x);times=int(sqrt(1+dy2)/sqrt(2*g*y),x=0.a);with(student);middlesum(sqrt(1+dy2)/sqrt(2*g*y),x=0.a,2000);Evalf(%);Time=0.57353,理论分析f(u)在u=0取得最小值, h(0)=0, h(a)=0, h(x) 光滑. f(0)=0.经过推导后得到方程利用 dsol
27、ve(y(x)*(diff(y(x),x)2+1)=c,y(x);求解后给出的结果太复杂,我们令,3.2 人体含铅量(第260页 19题) 铅是一种人体所必须的微量元素,但体内铅含量过多时就会引起铅中毒.铅是一种重金属元素,通过食物、饮料、空气进入人体,经过呼吸和消化系统后进入血液,再经过血液循环慢慢进入人体和骨头中.铅可以经过人体的排泄系统、通出出汗、剪头发、剪指甲排出体外. 我们根据铅在人体内的变化情况将人体分为血液、组织、骨头3个仓室,铅在这三个仓室中的转化关系如图所示.x1(t)表示t时刻血液中的含铅量,x2(t)表示t时刻组织中的含铅量,x3(t)表示t时刻骨头中的含铅量.假设在单位
28、时间内从环境经过消化、吸收系统进入血液的铅为L,从血液进入组织、骨头的铅分别为a31*x1(t)和a21*x1(t),从组织、骨头再进入血液的铅分别为a13*x2(t)、a12*x2(t),血液和骨头向外界排出的铅分别为a01*x1(t)和a02*x3(t). (1)根据前面的假设建立人体血液、组织和骨头中含铅量所满足的微分方程组,并求其通解.(2)取时间单位为天,对一个志愿者在一种环境中的生活情况进行测定得a21=0.011,a12=0.012, a31=0.0039,a13=0.000035,a01=0.021, a02=0.016,L=49.3. 假设该志愿者开始时体内的含铅量为0,求他
29、体内血液、组织和骨头中含铅量随时间变化的关系,画出这些解曲线的图,并求当x+时这些函数的极限.(3)假设该志愿者在此环境中生活了365天后搬到了一个无铅的环境中去(不再有外界的铅进入体内,即L=0),再讨论他体内血液、组织和骨头中含铅量随时间变化的关系,并画出0t1460时这些含铅量曲线的图形.,铅是一种重金属元素,在尘土、空中飘浮物、汽车尾气、油漆、一些彩釉陶瓷制品以及塑料所含的增塑剂中含量较高,主要通过口腔-消化道进入生物体,对其神经系统和血液系统有较强的毒性,会破坏其正常功能,其中以神经毒性为主。 众所周知,人的思想和行为受神经系统控制,一旦神经系统遭到破坏,其后果是可想而知的。儿童正处
30、于神经系统发育完善期,是铅污染的最大受害者。儿童好奇心强,喜欢用手去触摸各种物品,而且一般有“手-口接触(吮吸)”频繁的行为特点,这就给铅进入体内创造了条件。由于儿童正处于生长发育阶段,许多器官的功能都不完善,大脑的发育尚未完成,“血-脑屏障(一种保护性组织,有防止血液中物质随意进入脑细胞的功能)”不能阻止铅离子(2+)进入大脑。虽然,环境因素造成的铅危害一般不会达到工业性中毒的程度,但低量、长期的铅接触会造成铅在体内积累。虽然血铅含量超标并不等于铅中毒,但研究表明,过量的铅接触会造成神经元能量代谢降低,使儿童,尤其是学龄前儿童的神经发育指数达不到应有的标准,造成记忆力和思维能力下降,影响智力
31、的发展这种损害是不可逆的。,儿童体内的铅水平可分从血、骨、齿、尿、发检测到,儿童血铅水平分为5级。 级:血铅值低于10微克分升,身体处于相对安全状态; 级:血铅值为10微克19微克分升,属于轻度铅中毒。影响造血、神经传导和认知能力,儿童易出现头昏、烦躁、注意力涣散、多动、厌食、腹胀、轻度贫血; 级:血铅值达到20微克44微克分升,为中度铅中毒。引起缺钙、缺锌、缺铁、免疫力低下,运动不协调,视力和听力受损,学习困难、智力下降,生长发育迟缓,贫血、腹绞痛、食癖、反应迟钝等; 级:血铅值达到45微克69微克分升,就是重度铅中毒。可出现性格改变、易激怒、攻击性行为、运动失调、贫血、腹绞痛、高血压、心律
32、失常和痴呆等; 级:当血铅值大于70微克分升时,为极重度铅中毒,可导致脏器损害、铅性脑病、瘫痪、昏迷等。,铅中毒会引起神经、消化、循环系统紊乱,表现为贫血、腹痛、高血压等。血铅高组的儿童的总智商、操作智商、语言智商分别比低血铅组落后分、分和分,而每升血液中的铅浓度上升微克,儿童的身高将降低厘米。 低剂量的铅接触可以对人体的红细胞、肾脏、免疫系统、骨髓和中枢神经的功能产生不良影响,而所有这些影响发生前都可能没有明显的临床症状。铅中毒会影响婴幼儿最初站立、行走和说话的年龄,也可能引起孩子注意力涣散、记忆力减退、理解力降低及学习困难等。 铅通过各种方式进入人体后,可使人的智力下降,学习、工作成绩低落
33、;蓄积到一定程度时会使人出现精神障碍、噩梦、失眠、头痛等慢性中毒症状;严重者还可有乏力、食欲不振、恶心、腹胀、腹痛或腹泻等。,图,简化图,模型,模型x = Ax + b,数据,时间分为两段: 0,365和365,1460#定义常数a01:= 0.021; a02:=0.016; a12:=0.012; a13:=0.000035; a21:=0.011; a31:=0.0039; L:=49.3;#定义矩阵和向量A := matrix(3,3,-(a01+a21+a31), a12, a13, a21, -(a02+a12), 0, a31, 0, -a13);b :=matrix(3,1,L
34、, 0, 0);#定义方程equn1:=diff(x1(t),t)=-(a01+a21+a31)*x1(t)+a12*x2(t)+a13*x3(t)+L;equn2:=diff(x2(t),t)=a21*x1(t)-(a02+a12)*x2(t);equn3:=diff(x3(t),t)=a31*x1(t)-a13*x3(t);#解初始值问题dsolve(equn1,equn2,equn3,x1(0)=0,x2(0)=0,x3(0)=0,x1(t),x2(t),x3(t);,结果太复杂,利用迭加原理、特征值和特征向量法#非齐次方程的特解with(linalg): with(plots):xe
35、:= evalm(-(inverse(A),#得到解的表达式x1:=c1*P1,1*exp(lambda1*t)+c2*P1,2*exp(lambda2*t)+c3*P1,3*exp(lambda3*t)+xe1,1;x2 := c1*P2,1*exp(lambda1*t)+c2*P2,2*exp(lambda2*t)+c3*P2,3*exp(lambda3*t)+xe2,1;x3 := c1*P3,1*exp(lambda1*t)+c2*P3,2*exp(lambda2*t)+c3*P3,3*exp(lambda3*t)+xe3,1;#在解的表达式用t=0代入x10:=simplify(su
36、bs(t=0,x1);x20:=simplify(subs(t=0,x2);x30:=simplify(subs(t=0,x3);#求解常数solve(x10=0,x20=0,x30=0,c1,c2,c3);assign(%);,#作图plot(x1, x2, x3, t=0.365, color=red,blue,green, thickness=2); plot1:=%:,#重新显示解x1;x2;x3;#求极限limit(x1,t=infinity);limit(x2,t=infinity);limit(x3,t=infinity);#得到365天后解的表达式xx1:=cc1*P1,1*e
37、xp(lambda1*t)+cc2*P1,2*exp(lambda2*t)+cc3*P1,3*exp(lambda3*t);xx2 := cc1*P2,1*exp(lambda1*t)+cc2*P2,2*exp(lambda2*t)+cc3*P2,3*exp(lambda3*t);xx3 := cc1*P3,1*exp(lambda1*t)+cc2*P3,2*exp(lambda2*t)+cc3*P3,3*exp(lambda3*t);,# 计算解在t=365时的值 x1365:=simplify(subs(t=365,x1);x2365:=simplify(subs(t=365,x2);x3
38、365:=simplify(subs(t=365,x3);xx1365:=simplify(subs(t=365,xx1);xx2365:=simplify(subs(t=365,xx2);xx3365:=simplify(subs(t=365,xx3);#求解常数solve(xx1365=x1365,xx2365=x2365,xx3365=x3365, cc1,cc2,cc3);assign(%);#重新显示解xx1;xx2;xx3;#求极限limit(xx1,t=infinity);limit(xx2,t=infinity);limit(xx3,t=infinity);,#作图plot(xx1, xx2, xx3, t=365.1460, color=red,blue,green, thickness=2); plot2:=%:,#将两个图画在一个图上display(plot1,plot2);,