微分方程Sage快速入门

这个 Sage 快速入门教程是为MAA准备工作坊“Sage：对本科生使用开源数学软件”开发的（由NSF提供资金，到期日：0817071）。它是根据Creative Commons Attribution -ShareLiked 3.0许可证授权的 (CC BY-SA ）

求解微分方程是精确方法和数值方法的结合，因此是一个用计算机探索的好地方。我们已经在微积分教程中看到了一个例子，值得一看。

基本符号技巧

sage: y = function('y')(x)
sage: de = diff(y,x) + y -2
sage: h = desolve(de, y)

暂时忘记绘图，注意有三件事需要象征性地解微分方程：

  • 抽象函数 y ；

  • 一个微分方程，我们把它放在一条单独的直线上；

  • 实际 d 微分 e 方程 解决 命令（首字母缩写为粗体 desolve ）

sage: show(expand(h))
\[c e^{左（-x右）}+2\]

由于我们没有指定任何初始条件，Sage（来自Maxima）输入了一个参数。如果我们想进入初始状态，我们使用 ics （用于 i 初始 c 条件s ). 例如，我们设置 ics=[0,3] 指定何时 x=0 ， y=3 .

sage: h = desolve(de, y, ics=[0,3]); h
(2*e^x + 1)*e^(-x)

当然，我们已经注意到，我们可以用斜率场来绘制所有这些。

sage: y = var('y') # Needed so we can plot
sage: Plot1=plot_slope_field(2-y,(x,0,3),(y,0,5))
sage: Plot2=plot(h,x,0,3)
sage: Plot1+Plot2
Graphics object consisting of 2 graphics primitives

注解

关于符号函数与符号变量：

  • 如果你想 y 又是一个抽象函数而不是一个变量，你必须分开来做。微分方程需要 function 但策划需要 var .

  • 另一个选择是 z 是垂直轴变量的名称。

  • 不管怎样，我们都要付出，因为我们对待这些事情的共同点是 y 作为一个变量和一个函数，这要用计算机来完成要困难得多。

Sage中还有许多其他的微分方程工具。我们不能在快速入门中涵盖所有的变体，但是文档对符号求解器很有用。

sage: desolvers?

例如，Maxima可以做系统，也可以使用Laplace变换，为了便于使用，我们将这些封装的版本包括在内。

在所有的微分方程求解程序中，注意语法是很重要的！在下面的示例中，我们将微分方程放在命令的主体中，并且必须指定 f 是 d 依附 var 可行的 (dvar )，并给出初始条件 \(f(0)=1\) 和 \(f'(0)=2\) ，它给出了示例中的最后一个列表。

sage: f=function('f')(x)
sage: desolve_laplace(diff(f,x,2) == 2*diff(f,x)-f, dvar = f, ics = [0,1,2])
x*e^x + e^x

sage: g(x)=x*e^x+e^x
sage: derivative(g,x,2)-2*derivative(g,x)+g
x |--> 0

数值和幂级数方法

也有数值方法。

例如，上面的选项之一是 desolve_rk4 . 这是一种四阶龙格库塔方法，并返回适当的（数值）输出。在这里，我们 must 给出因变量 and 初始条件。

sage: y = function('y')(x)
sage: de = diff(y,x) + y -2
sage: h = desolve_rk4(de, y, step=.05, ics=[0,3])

将其与原始的象征性解决方案进行比较会很有趣。我们使用 points 高级打印教程中的命令。

sage: h1 = desolve(de, y, ics=[0,3])
sage: plot(h1,(x,0,5),color='red')+points(h)
Graphics object consisting of 2 graphics primitives

从这里开始，数字例程的主要用途是教学性的。

对于更高级的数值程序，我们主要使用GNU科学库。使用它稍微复杂一点，但提供了丰富的选择。

sage: ode_solver?

我们甚至可以做幂级数解。为了做到这一点，我们必须首先定义一个特殊的 幂级数环 包括精度。

sage: R.<t> = PowerSeriesRing(QQ, default_prec=10)
sage: a = -1 + 0*t
sage: b = 2 + 0*t
sage: h=a.solve_linear_de(b=b,f0=3,prec=10)

这个幂级数解决方案暂时还不错！

sage: h = h.polynomial()
sage: plot(h,-2,5)+plot(2+e^-x,(x,-2,5),color='red',linestyle=':',thickness=3)
Graphics object consisting of 2 graphics primitives

这只是一个介绍；在其他地方有很多使用Sage的微分方程的资源，包括davidjoyner的一本书，他为Sage编写了很多包装Maxima的原始代码。