1 / 20第四章 微分方程§4.1 常微分方程4.1.1 常微分方程的解析解1. 函数dsolve 在微分方程中的应用在Maple 中,这是一个用途最广的函数——称为通用函数吧,几乎可以求解所有的微 分方程和方程组,既能求解解析解,也能求解数值解,本节只介绍其求微分方程的解析解中的作用:dsolve (ODE);dsolve (ODE,y(x),extra_args);其中,ODE(Ordinary Differential Equation)是一个常微分方程; y(x)为未函数,求解时这个参数可以省略;第三个参数extra_args 是一个可选的参数,主要用来设置最后解析解的形式或求解过程中一些积分的设置,它的选值很广,这里仅举几个参数。
(1) explicit: 求出显式解; (2) implicit: 解可以是隐式;(3) useInt: 运算中用“Int ”函数代替“int ”函数,可加快运算速度; (4) parametric: 将最后的解析解表达成另外一个自变量的形式。
这些参数的位置很灵活,可以放在除第一个参数位置外的任何位置,并且它们的组合 也很灵活,可以单独作用,也可几何参数合用,只要在中间用逗号隔开,而且参数并不一定需要写在一起,也可以分开。
> eq:='eq': eq:=diff(y(x),x)*(1+y(x)^2)+cos(x)=0; 可以两端都不是零:= eq = + ⎛⎝ ⎫⎭⎪⎪∂()y x () + 1()y x 2()cos x 0 > sol1:=dsolve(eq,explicit); 给出显式解sol1()y x =:= 12- ()- - + 12()sin x 12_C12 - + + 3418()cos 2x 72()sin x _C136_C12()/234()- - + 12()sin x 12_C12 - + + 3418()cos 2x 72()sin x _C136_C12()/13,其中“_C 1”表示第一个任意常数。
方程的解是很恐怖的解,这里仅给出了一个解,另2 / 20外还有两个更长的解,读者可以在Maple 下执行上面求解过程观察到另外两个解的全貌。
这是由于将解转换成显函数造成的,假如我们将参数进行改善:> sol2:=dsolve(eq,implicit,y(x)); 给出隐式解,式中的y(x)可省略:= sol2 = + + + ()sin x ()y x 1()y x 3_C10再加上一个参数“useInt ”,可以明显感到运算速度非常快,因此,它在求解过程中积分比较复杂时很有用,同时还能使解过程、解结果给出较多的信息: > sol3:=dsolve(eq,implicit,y(x),useInt);:= sol3 = - + d ⎛⎠⎜()cos x x d ⎛⎠⎜⎜()y x - - 1_a 2_a _C10其中“_a ”为积分变量.即解为.0)1(cos 2=+---⎰⎰C dt t xdx ya最后加入参数“parametric ”,可以知道经过一段时段运算后的结果: > sol4:=dsolve(eq,implicit,y(x),useInt,parametric);:= sol4我们惊讶地发现函数没有给出任何结果,这是因为解太复杂了,函数找不到用参数表示的方法。
下面我们用一个比较简单的例子来说明设置参数以后的结果,大家容易从结果中看出表示的方法:>dsolve(diff(y(x),x)=-x/y(x), parametric); 圆曲线上切线的斜率, = ()y x - + x 2_C1 = ()y x -- + x 2_C1此为圆的上半圆与下半圆曲线表示式.> dsolve(diff(y(x),x)=-x/y(x),implicit, parametric); 参数式⎡⎣⎢⎢⎤⎦⎥⎥, = ()x _T -_T _C1 + 1_T 2 = ()y _T _C1 + 1_T 2此为圆曲线的参数式,但并不是常用的t a y t a x sin ,cos ==参数式格式. 2. 用函数odetest 检验常微分方程的解odetest(sol,ODE); ——y(x)可省略 odetest(sol,ODE,y(x));——y(x)最好加上 odetest(solsys,sysODE);——用于方程组以返回值为“0”给出解为真。
> with(DEtools):3 / 20> odetest(sol1[1],eq,y(x)); sol1[1]是方程的解> odetest(sol1[2],eq,y(x)); sol1[2]是方程的解> odetest(sol1[3],eq,y(x)); sol1[3]是方程的解> odetest(sol2,eq,y(x)); sol2是方程的解> odetest(sol3,eq,y(x)); sol3是方程的解> odetest(sol4,eq,y(x)); sol4不能代入检验Error, (in odetest) expecting the second argument to be an ODE or a set or list of ODEs. Received: y(x)下面验证一个函数是否前面所给方程的解:> y(x)=x^2; 验证y(x)=x^2是否方程的解,这里的y(x)不能赋给= ()y x x 2> odetest(%,eq); 所给函数y(x)=x^2使方程左端不为0,故不是方程的解+ + 2x 2x 5()cos x但它是下列方程的解:> eq:='eq': eq:=diff(y(x),x)=2*x;:= eq = ∂()y x 2x> y(x)=x^2;= ()y x x 2> odetest(%,eq);3.用Deplot 函数来显示微分方程的解的图像DEplot(deqns,vars,trange,eqns); Deplot(deqns,vars,trange,inits, eqns);其中 deqns :是待求的微分方程或微分方程组,本节中就是指微分方程;vars :是一个非独立的变量或变量的集合,如y(x);trange :设定变量的取值范围;eqns 是一个等式形式的条件,如:linecolor=blue(设定图像中线的颜色为蓝的);4 / 20stepsize=0.05(变量的取值步长);inits :是微分方程的初值条件,如果有多个初值条件,用集合表示,并且每个初值条件都必须用方括号起来,如果不给初始值,函数将等间隔的取最后解中的常数为某值。
请看下面的例子:> with(DEtools):> eq:='eq': eq:=diff(y(x),x)*(1+y(x)^2)+cos(x)=0;:= eq = + ⎛⎝ ⎫⎭⎪⎪∂∂x ()y x () + 1()y x 2()cos x 0 > DEplot(eq,y(x),x=-10..10,y=0..3,linecolor=blue,stepsize=0.05);这个例子中没有为微分方程设定初始值,所以函数用一些箭头来表示解的趋势,假如我们设方程的初始值是y(0)=2。
图像将变为:> DEplot(eq,y(x),x=-10..10,[[y(0)=2]],y=0..3,linecolor=blue,stepsize =0.05); [[y(0)=2]]这个初始条件也可放在“y=0..3”的后面图形中的粗线就是满足初始条件的解的图形,假如我们对这个图形是不是刚才求得的解的图形有疑问的话,我们可以用另外一个作图函数contourplot 来验证一下,这个函数在plots 程序包里: > with(plots):5 / 20> f:=sin(x)+1/3*y^3+y:> contourplot(f,x=-10..10,y=-3..3,coloring=[red,blue]);显然,过y(0)=2这点的图形跟上面图形中篮粗线对应的图形是一样的,说明我们的求解没有错。
如果将这个图的范围适当调整,尽量得到过y(0)=2的那条曲线,显示效果会更好些。
实际上,只要将作图中y 的范围改为-2.55~2.55即可达到目的.4.用各种特定函数求解各种微分方程——P176odeadvisor(ODE); 在Detools 程序包中这是一个判断微分方程类型的函数,可用来判分离变量、一阶线性、全微分方程还是普通方程。
最常见的返回参数和中文含义如下:separable 可分离变量方程,用separablesol (lode,v)求解,其中v 为未知函数 linear 线性方程, 用linearsol (lode,v)homogeneous 齐次型方程,用homogeneous (lode,v) bernoulli 贝努里 exact 全微分方程 Abel 阿贝尔方程Quadrature 积分形式的微分方程(如dy=xdx ),是可离变量方程中最简单的(1)可分离变量的方程.)1( Ex.1y y x y x xy '+='-解可分离变量方程解:(Maple-lx4-解法参考) > eq:='eq':> eq:=x*y(x)*(1-x*diff(y(x),x))=x+y(x)*diff(y(x),x);:= eq = x ()y x ⎛⎝ ⎫⎭⎪⎪ - 1x ⎛⎝ ⎫⎭⎪⎪∂∂x ()y x + x ()y x ⎛⎝ ⎫⎭⎪⎪∂∂x ()y x > with(DEtools): 调用程序包6 / 20> odeadvisor(eq); 判断方程类型[]_separable 方程为可分离变量的方程> separablesol(eq,y(x)); 用程序包中的分离变量法求解{} = ()y x + ()LambertW + x 21_C1e()-11法2: 仍用dsolve 求解:> dsolve(eq,explicit); 用通用函数求解.接上方程,结果与上法相同= ()y x + ()LambertW + x 21_C1e()-11法3: 将方程改为C dx x f dy y g y g x f dx dy +==⎰⎰)()(,)()(再利用求隐式解: > diffy:='diffy': 初始化变量,这步一要定做 > diffy=solve(eq,diff(y(x),x)); 解出y 的导数=diffy x ()- ()y x 1()y x ()+ x 21 > f:='f': g:='g':> f:=x->x/(x^2+1);g:=y->y/(y-1); 构造函数f(x),g(y):= f →x x+ x 21:= g →y y这里的y 视为单变量,不要写成y(x)> RHS:=int(g(y),y);LHS:=int(f(x),x); 赋表达式分别得左、右两端:= RHS + y ()ln - y 1:= LHS 12()ln + x 21> LHS=RHS+C; 给出隐式解= 12()ln + x 21 + + y ()ln - y 1C > solve(%,y); 求出显式解+ ()LambertW e()- - /12()ln + x 211C 17 / 20(2)齐次型的方程.tan Ex.2的通解求齐次型方程xyx y y +=' 解 (Maple-lx4-1)法1: (按齐次方程专用求解函数求解) > eq:='eq';:= eq eq> eq:=diff(y(x),x)=y(x)/x+tan(y(x)/x);:= eq = ∂∂x ()y x + ()y x x ⎛⎝ ⎫⎭⎪⎪tan ()y x x > with(DEtools): 调程序包 > odeadvisor(eq); 判断类型[],,[],_homogeneous class A [],_1st_order _with_linear_symmetries _dAlembert齐次型方程的A 种,是一阶、线性对称方程> genhomosol(eq); 用齐次方程的求解函数求解,得以下两个通解{}, = ()y x ⎛⎝⎫⎭⎪⎪⎪arctan -()- + 1_C1x 2_C1x - + 1_C1x 2x = ()y x -⎛⎝ ⎫⎭⎪⎪⎪arctan -()- + 1_C1x 2_C1x - + 1_C1x 2x 法2: 用通用函数求解:> dsolve(eq,explicit);, = ()y x ⎛⎝ ⎫⎭⎪⎪⎪arctan -()- + 1_C1x 2_C1x - + 1_C1x 2x = ()y x -⎛⎝ ⎫⎭⎪⎪⎪arctan -()- + 1_C1x 2_C1x - + 1_C1x 2x法3: 按解齐次型方程的求解方法求解:> y:='y':x:='x':u:='u': 如果不是新打开的Maple 工作窗口,最好要这一步 > y:=u*x:> one:=expand(eq); 试将方程中“y/x”直接换成u 后给出新方程:= one = + ⎛⎝ ⎫⎭⎪⎪∂∂x ()u x ()x x ()u x ⎛⎝ ⎫⎭⎪⎪∂∂x ()x x + ()u x ()x x x ⎛⎝ ⎫⎭⎪⎪tan ()u x ()x x x > two:=subs(diff(x(x),x)=1,x(x)=x,one); 简化方程:= two = + ⎛⎝ ⎫⎭⎪⎪∂∂x ()u x x ()u x + ()u x ()tan ()u x8 / 20> diffu:=solve(two,diff(u(x),x)); 化为可分离变量的方程:=diffu ()tan ()u x x> left:=int(1/tan(u),u);right:=int(1/x,x); 按可分离变量方法求解,这里的u 不能写成u(x),因为要把它当单变量进行积分:= left - ()ln ()tan u 12()ln + 1()tan u 2:= right ()ln x> y:='y': 初始化变量y,此前它是”ux” > subs(u=y/x,left)=right+_C1;= - ⎛⎝ ⎫⎭⎪⎪ln ⎛⎝ ⎫⎭⎪⎪tan y x 12⎛⎝ ⎫⎭⎪⎪⎪ln + 1⎛⎝ ⎫⎭⎪⎪tan y x 2+ ()ln x _C1 > solve(%,y);,⎛⎝⎫⎭⎪⎪⎪⎪arctan -()- + 1e()2_C1x 2e()2_C1x - + 1e()2_C1x 2x -⎛⎝⎫⎭⎪⎪⎪⎪arctan -()- + 1e()2_C1x 2e()2_C1x - + 1e()2_C1x 2x 这一结果与前的解仅在于任意常数的写法上不同而已.法4: 本人的解法——化为分离变量方程的求解函数进行求解:> y(x):='y(x)':x:='x':u(x):='u(x)':> y(x):=u(x)*x: 注意函数关系,如果用y=ux 或(x)=u*x,得不到下面的形式> eq1:=subs(y(x)/x=u(x),eq); 替换完毕,再对方程化简:= eq1 = + ⎛⎝ ⎫⎭⎪⎪∂∂x ()u x x ()u x + ()u x ()tan ()u x上式两端可不再化简,解过程用时相近> separablesol(eq1,u(x)); 用分离变量法函数求解方程,得两个同样的通解{}, = ()u x ⎛⎝ ⎫⎭⎪⎪⎪arctan () - 1_C1x 2_C1x - + 1_C1x 2 = ()u x -⎛⎝ ⎫⎭⎪⎪⎪arctan () - 1_C1x 2_C1x - + 1_C1x 2 > a:=subs(u(x)=y/x,%); 换回原变量,如果不将两个通解赋值给集合a,就不能用“solve ”对集合中的两个非方程组的函数关于同一变量求解9 / 20:= a {}, = y ⎛⎝ ⎫⎭⎪⎪⎪arctan -()- + 1_C1x 2_C1x - + 1_C1x 2 = y -⎛⎝ ⎫⎭⎪⎪⎪arctan -()- + 1_C1x 2_C1x - + 1_C1x 2 > y1:=solve(op(1,a),y);y2:=solve(op(2,a),y); 给出两个显式解:= y1⎛⎝ ⎫⎭⎪⎪⎪arctan -()- + 1_C1x 2_C1x - + 1_C1x 2x := y2-⎛⎝⎫⎭⎪⎪⎪arctan -()- + 1_C1x 2_C1x - + 1_C1x 2x 上面用了提出集合元素的函数“op ”,参第一章“1.1.3 集合”. 也可用分离方程中变量的函数isolate,执行“isolate(op(1,a),y);”得到第一个通解中的y 的表达式.请试之.Ex.2’ 解齐次型(C )的微分方程.31-++-='y x y x y 答案: = - + - - y 24y 2x y 2x x 2C解 (Maple-lx4-2) .,)(222111解读者可比照固定方法求有固定的解法对准齐次型方程c y b x a c y b x a f dx dy++++=也可以用通用求解函数dsolve()求解.这里选用求齐次型方程的求解函数进行求解: > eq:='eq':> eq:=diff(y(x),x)=(x-y(x)+1)/(x+y(x)-2);:= eq = ∂∂x ()y x - + x ()y x 1 + - x ()y x 2> with(DEtools):> a:=genhomosol(eq); 用齐次型方程的求解函数求解,解用集合表示,其解需要按元素取出使用,为此,将解赋给集合“a ”:= a {} = ()y x - 31 + () - 2x 1_C1 + 2() - 2x 12_C121由于x 与y 的隐式关系很简单,而这里的显示结果即较为复杂,为此,如以下化简:> subs(_C1=C,a[1]):subs(y(x)=z,a[1]):subs(_C1=C,%); 替换= z - 3212 + () - 2x 1C + 2() - 2x 12C 21C10 / 20> map(simplify,%); 对上式两端化简。