编程

加载和附加Sage文件

接下来,我们将说明如何将以单独文件编写的程序加载到Sage中。创建名为的文件 example.sage 内容如下:

print("Hello World")
print(2^3)

您可以读入并执行 example.sage 文件中使用 load 指挥部。

sage: load("example.sage")
Hello World
8

方法将Sage文件附加到正在运行的会话 attach 命令:

sage: attach("example.sage")
Hello World
8

现在如果你改变了 example.sage 并在Sage中输入一个空行(即,点击 return ),则内容 example.sage 将自动重新加载到Sage中。

特别是, attach 每当文件更改时自动重新加载文件,这在调试代码时很方便,而 load 一个文件只加载一次。

当Sage加载时 example.sage 它将其转换为Python,然后由Python解释器执行。这种转换是最小的;它主要涉及在 Integer() 中的浮点文字 RealNumber() ,替换 ^ 是由 ** 的,并替换,例如, R.2 通过 R.gen(2) 。的转换版本 example.sage 包含在与相同的目录中 example.sage 并被称为 example.sage.py 。该文件包含以下代码:

print("Hello World")
print(Integer(2)**Integer(3))

整型文字被包装起来,并且 ^ 被替换为 ** 。(在Python中 ^ 意思是“异或”和 ** 意思是“指数运算”。)

(此准备工作在 sage/misc/interpreter.py 。)

您可以将多行缩进代码粘贴到Sage中,只要有换行符来创建新块(这在文件中不是必需的)。然而,将这样的代码输入到Sage中的最佳方式是将其保存到文件中并使用 attach ,如上所述。

创建编译后的代码

速度在数学计算中是至关重要的。尽管Python是一种方便的高级语言,但如果某些计算是使用编译语言中的静态类型实现的,那么它们可能会比在Python语言中快几个数量级。如果Sage完全用Python语言编写,那么它的某些方面可能会太慢。为了解决这一问题,Sage支持名为Cython的已编译版本的Python ([Cyt][Pyr]) 。Cython同时类似于Python和C语言。大多数Python构造,包括列表理解、条件表达式、如下代码 += 您还可以导入在其他Python模块中编写的代码。此外,您可以声明任意的C变量,并且可以直接进行任意的C库调用。结果代码被转换为C语言,并使用C编译器进行编译。

为了生成您自己的已编译的Sage代码,请给该文件一个 .spyx 扩展(而不是 .sage )。如果您使用的是命令行界面,则可以像使用解释代码一样附加和加载编译后的代码(目前,笔记本界面不支持附加和加载Cython代码)。实际的编译是在“幕后”完成的,您不必进行任何显式操作。编译后的共享对象库存储在 $HOME/.sage/temp/hostname/pid/spyx 。退出Sage时,这些文件将被删除。

不对SPYX文件应用SAGE准备,例如, 1/3 将导致SPYX文件中的0而不是有理数 \(1/3\) 。如果 foo 是Sage库中的一个函数,用于从SPYX文件导入使用它 sage.all 和使用 sage.all.foo

import sage.all
def foo(n):
    return sage.all.factorial(n)

在单独的文件中访问C函数

访问在单独的*.c文件中定义的C函数也很容易。这里有一个例子。创建文件 test.ctest.spyx 在包含以下内容的同一目录中:

纯C代码: test.c

int add_one(int n) {
  return n + 1;
}

Cython代码: test.spyx

cdef extern from "test.c":
    int add_one(int n)

def test(n):
    return add_one(n)

然后,下面的方法就会奏效:

sage: attach("test.spyx")
Compiling (...)/test.spyx...
sage: test(10)
11

如果一个额外的库 foo 需要编译从Cython文件生成的C代码,添加下面这行 clib foo 添加到Cython源代码中。类似地,另一个C文件 bar 可以与声明一起包含在编译中 cfile bar

独立的Python/Sage脚本

以下独立的Sage脚本将整数、多项式等作为因子:

#!/usr/bin/env sage

import sys
from sage.all import *

if len(sys.argv) != 2:
    print("Usage: %s <n>" % sys.argv[0])
    print("Outputs the prime factorization of n.")
    sys.exit(1)

print(factor(sage_eval(sys.argv[1])))

为了使用此脚本,您的 SAGE_ROOT 一定是挡住了你的去路。如果调用上述脚本 factor ,下面是一个用法示例:

$ ./factor 2006
2 * 17 * 59

数据类型

Sage中的每个对象都有一个定义良好的类型。Python具有广泛的基本内置类型,Sage库还添加了更多类型。一些内置的Python类型包括字符串、列表、元组、整型和浮点型,如下所示:

sage: s = "sage"; type(s)
<... 'str'>
sage: s = 'sage'; type(s)      # you can use either single or double quotes
<... 'str'>
sage: s = [1,2,3,4]; type(s)
<... 'list'>
sage: s = (1,2,3,4); type(s)
<... 'tuple'>
sage: s = int(2006); type(s)
<... 'int'>
sage: s = float(2006); type(s)
<... 'float'>

除此之外,Sage还添加了许多其他类型的内容。例如,向量空间:

sage: V = VectorSpace(QQ, 1000000); V
Vector space of dimension 1000000 over Rational Field
sage: type(V)
<class 'sage.modules.free_module.FreeModule_ambient_field_with_category'>

只能调用某些函数 V 。在其他数学软件系统中,这些函数将使用“函数”符号来调用 foo(V,...) 。在Sage中,某些函数附加到类型(或类) V ,并且使用类似于Java或C++的面向对象的语法来调用,例如, V.foo(...) 。这有助于防止全局命名空间被数万个函数污染,这意味着许多具有不同行为的不同函数可以被命名为foo,而不必使用参数(或case语句)的类型检查来决定调用哪个函数。此外,如果您重复使用某个函数的名称,该函数仍然可用(例如,如果您调用了某个函数 zeta ,然后想要计算Riemann-Zeta函数的值为0.5,您仍然可以输入 s=.5; s.zeta() )。

sage: zeta = -1
sage: s=.5; s.zeta()
-1.46035450880959

在一些非常常见的情况下,为了方便起见,还支持常用的函数表示法,因为使用面向对象的表示法可能会使数学表达式看起来令人困惑。这里有一些例子。

sage: n = 2; n.sqrt()
sqrt(2)
sage: sqrt(2)
sqrt(2)
sage: V = VectorSpace(QQ,2)
sage: V.basis()
    [
    (1, 0),
    (0, 1)
    ]
sage: basis(V)
    [
    (1, 0),
    (0, 1)
    ]
sage: M = MatrixSpace(GF(7), 2); M
Full MatrixSpace of 2 by 2 dense matrices over Finite Field of size 7
sage: A = M([1,2,3,4]); A
[1 2]
[3 4]
sage: A.charpoly('x')
x^2 + 2*x + 5
sage: charpoly(A, 'x')
x^2 + 2*x + 5

列出的所有成员函数 \(A\) ,使用制表符完成。只需输入即可 A. ,然后键入 [tab] 键盘上的键,如中所述 反向搜索和制表符完成

列表、元组和序列

列表数据类型存储任意类型的元素。与C、C++等类似(但与大多数标准计算机代数系统不同),列表的元素从 \(0\)

sage: v = [2, 3, 5, 'x', SymmetricGroup(3)]; v
[2, 3, 5, 'x', Symmetric group of order 3! as a permutation group]
sage: type(v)
<... 'list'>
sage: v[0]
2
sage: v[2]
5

(在索引到列表时,如果索引不是Pythonint!)Sage整数(或有理数,或带有 __index__ 方法)将工作得很好。

sage: v = [1,2,3]
sage: v[2]
3
sage: n = 2      # Sage Integer
sage: v[n]       # Perfectly OK!
3
sage: v[int(n)]  # Also OK.
3

这个 range 函数创建一组Python int(不是Sage整数):

sage: list(range(1, 15))
[1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14]

在使用列表理解构造列表时,这很有用:

sage: L = [factor(n) for n in range(1, 15)]
sage: L
[1, 2, 3, 2^2, 5, 2 * 3, 7, 2^3, 3^2, 2 * 5, 11, 2^2 * 3, 13, 2 * 7]
sage: L[12]
13
sage: type(L[12])
<class 'sage.structure.factorization_integer.IntegerFactorization'>
sage: [factor(n) for n in range(1, 15) if is_odd(n)]
[1, 3, 5, 7, 3^2, 11, 13]

有关如何使用列表理解创建列表的更多信息,请参见 [PyT].

列表切片是一个很棒的功能。如果 L 是一份名单,那么 L[m:n] 返回的子列表 L 通过从 \(m^{th}\) 元素,并在 \((n-1)^{st}\) 元素,如下所示。

sage: L = [factor(n) for n in range(1, 20)]
sage: L[4:9]
[5, 2 * 3, 7, 2^3, 3^2]
sage: L[:4]
[1, 2, 3, 2^2]
sage: L[14:4]
[]
sage: L[14:]
[3 * 5, 2^4, 17, 2 * 3^2, 19]

元组类似于列表,除了它们是不可变的,这意味着一旦它们被创建,它们就不能被改变。

sage: v = (1,2,3,4); v
(1, 2, 3, 4)
sage: type(v)
<... 'tuple'>
sage: v[1] = 5
Traceback (most recent call last):
...
TypeError: 'tuple' object does not support item assignment

序列是第三种面向列表的Sage类型。与列表和元组不同,Sequence不是内置的Python类型。默认情况下,序列是可变的,但使用 Sequence 类方法 set_immutable ,可以将其设置为不可变,如下面的示例所示。序列的所有元素都有一个共同的父元素,称为序列全域。

sage: v = Sequence([1,2,3,4/5])
sage: v
[1, 2, 3, 4/5]
sage: type(v)
<class 'sage.structure.sequence.Sequence_generic'>
sage: type(v[1])
<class 'sage.rings.rational.Rational'>
sage: v.universe()
Rational Field
sage: v.is_immutable()
False
sage: v.set_immutable()
sage: v[0] = 3
Traceback (most recent call last):
...
ValueError: object is immutable; please change a copy instead.

序列派生自列表,可以在任何可以使用列表的地方使用:

sage: v = Sequence([1,2,3,4/5])
sage: isinstance(v, list)
True
sage: list(v)
[1, 2, 3, 4/5]
sage: type(list(v))
<... 'list'>

作为另一个例子,向量空间的基是不变的序列,因为重要的是不要改变它们。

sage: V = QQ^3; B = V.basis(); B
[
(1, 0, 0),
(0, 1, 0),
(0, 0, 1)
]
sage: type(B)
<class 'sage.structure.sequence.Sequence_generic'>
sage: B[0] = B[1]
Traceback (most recent call last):
...
ValueError: object is immutable; please change a copy instead.
sage: B.universe()
Vector space of dimension 3 over Rational Field

辞典

字典(有时也称为关联数组)是从“Hasable”对象(例如,字符串、数字和此类对象的元组;有关详细信息,请参阅http://docs.python.org/tut/node7.html和http://docs.python.org/lib/typesmapping.html文档)到任意对象的映射。

sage: d = {1:5, 'sage':17, ZZ:GF(7)}
sage: type(d)
<... 'dict'>
sage: list(d.keys())
[1, 'sage', Integer Ring]
sage: d['sage']
17
sage: d[ZZ]
Finite Field of size 7
sage: d[1]
5

第三个关键字说明词典的索引可能很复杂,例如整数环。

您可以将上面的字典转换为包含相同数据的列表:

sage: list(d.items())
[(1, 5), ('sage', 17), (Integer Ring, Finite Field of size 7)]

一个常见的习语是遍历词典中的对:

sage: d = {2:4, 3:9, 4:16}
sage: [a*b for a, b in d.items()]
[8, 27, 64]

正如最后一条输出所示,词典是无序的。

Python有一个内置的set类型。它提供的主要功能是非常快速地查找元素是否在集合中,以及标准集合论操作。

sage: X = set([1,19,'a']);   Y = set([1,1,1, 2/3])
sage: X   # random sort order
{1, 19, 'a'}
sage: X == set(['a', 1, 1, 19])
True
sage: Y
{2/3, 1}
sage: 'a' in X
True
sage: 'a' in Y
False
sage: X.intersection(Y)
{1}

Sage也有自己的Set类型,该类型(在某些情况下)是使用内置的Python Set类型实现的,但具有一些额外的Sage相关功能。使用创建Sage Set Set(...) 。例如,

sage: X = Set([1,19,'a']);   Y = Set([1,1,1, 2/3])
sage: X   # random sort order
{'a', 1, 19}
sage: X == Set(['a', 1, 1, 19])
True
sage: Y
{1, 2/3}
sage: X.intersection(Y)
{1}
sage: print(latex(Y))
\left\{1, \frac{2}{3}\right\}
sage: Set(ZZ)
Set of elements of Integer Ring

迭代器

迭代器是在数学应用程序中特别有用的新近添加到Python中的。以下是几个示例;请参见 [PyT] 了解更多详细信息。我们在非负整数的平方上做一个迭代器,直到 \(10000000\)

sage: v = (n^2 for n in range(10000000))
sage: next(v)
0
sage: next(v)
1
sage: next(v)
4

我们在形式的素数上创建一个迭代 \(4p+1\) 使用 \(p\) 也是质数,看看前几个值。

sage: w = (4*p + 1 for p in Primes() if is_prime(4*p+1))
sage: w         # in the next line, 0xb0853d6c is a random 0x number
<generator object at 0xb0853d6c>
sage: next(w)
13
sage: next(w)
29
sage: next(w)
53

某些环(例如有限域和整数)具有与其相关联的迭代器:

sage: [x for x in GF(7)]
[0, 1, 2, 3, 4, 5, 6]
sage: W = ((x,y) for x in ZZ for y in ZZ)
sage: next(W)
(0, 0)
sage: next(W)
(0, 1)
sage: next(W)
(0, -1)

循环、函数、控制语句和比较

我们已经看到了一些常见用法的例子 for 循环。在Python中,一个 for 循环具有缩进结构,例如

>>> for i in range(5):
...     print(i)
...
0
1
2
3
4

注意for语句末尾的冒号(在Gap或Maple中没有“do”或“od”),以及循环“Body”之前的缩进,即 print(i) 。这种缩进很重要。在Sage中,当您点击时,缩进会自动放入 enter 在“:”之后,如下所示。

sage: for i in range(5):
....:     print(i)  # now hit enter twice
....:
0
1
2
3
4

这个符号 = 用于赋值。这个符号 == 用于检查等价性:

sage: for i in range(15):
....:     if gcd(i,15) == 1:
....:         print(i)
....:
1
2
4
7
8
11
13
14

请记住,缩进如何决定 iffor ,以及 while 声明:

sage: def legendre(a,p):
....:     is_sqr_modp=-1
....:     for i in range(p):
....:         if a % p == i^2 % p:
....:             is_sqr_modp=1
....:     return is_sqr_modp

sage: legendre(2,7)
1
sage: legendre(3,7)
-1

当然,这不是Legendre符号的有效实现!它旨在说明Python/Sage编程的各个方面。Sage附带的函数{kronecker}通过对pari的C库调用高效地计算Legendre符号。

最后,我们注意到比较,例如 ==!=<=>=>< ,如果可能,Better Numbers将自动将这两个数字转换为同一类型:

sage: 2 < 3.1; 3.1 <= 1
True
False
sage: 2/3 < 3/2;   3/2 < 3/1
True
True

使用bool表示符号不等:

sage: x < x + 1
x < x + 1
sage: bool(x < x + 1)
True

在比较Sage中不同类型的对象时,在大多数情况下,Sage会尝试找到两个对象对共同父对象的规范强制(请参见 父母、皈依和胁迫 以了解更多详细信息)。如果成功,则在强制对象之间执行比较;如果不成功,则认为这些对象不相等。要测试两个变量是否引用相同的对象,请使用 is 。正如我们在本例中看到的,Pythonint 1 是唯一的,但Sage整数 1 不是:

sage: 1 is 2/2
False
sage: 1 is 1
False
sage: 1 == 2/2
True

In the following two lines, the first equality is False because there is no canonical morphism \(\QQ\to \GF{5}\), hence no canonical way to compare the \(1\) in \(\GF{5}\) to the \(1 \in \QQ\). In contrast, there is a canonical map \(\ZZ \to \GF{5}\), hence the second comparison is True. Note also that the order doesn't matter.

sage: GF(5)(1) == QQ(1); QQ(1) == GF(5)(1)
False
False
sage: GF(5)(1) == ZZ(1); ZZ(1) == GF(5)(1)
True
True
sage: ZZ(1) == QQ(1)
True

WARNING: Comparison in Sage is more restrictive than in Magma, which declares the \(1 \in \GF{5}\) equal to \(1 \in \QQ\).

sage: magma('GF(5)!1 eq Rationals()!1')            # optional - magma
true

仿形

“过早优化是万恶之源。”--唐纳德·克努斯

有时,检查代码中的瓶颈以了解哪些部分占用了最多的计算时间是很有用的;这可以很好地确定要优化哪些部分。因此,Python和Sage提供了几个分析--这个过程被称为--选项。

最简单的用法是 prun 命令在交互Shell中。它返回一个摘要,描述哪些函数花费了多少计算时间。例如,要分析有限域上的矩阵乘法(当前最慢!-从1.0版开始),请执行以下操作:

sage: k,a = GF(2**8, 'a').objgen()
sage: A = Matrix(k,10,10,[k.random_element() for _ in range(10*10)])
sage: %prun B = A*A
       32893 function calls in 1.100 CPU seconds

Ordered by: internal time

ncalls tottime percall cumtime percall filename:lineno(function)
 12127  0.160   0.000   0.160  0.000 :0(isinstance)
  2000  0.150   0.000   0.280  0.000 matrix.py:2235(__getitem__)
  1000  0.120   0.000   0.370  0.000 finite_field_element.py:392(__mul__)
  1903  0.120   0.000   0.200  0.000 finite_field_element.py:47(__init__)
  1900  0.090   0.000   0.220  0.000 finite_field_element.py:376(__compat)
   900  0.080   0.000   0.260  0.000 finite_field_element.py:380(__add__)
     1  0.070   0.070   1.100  1.100 matrix.py:864(__mul__)
  2105  0.070   0.000   0.070  0.000 matrix.py:282(ncols)
  ...

这里 ncalls 是呼叫数, tottime 是在给定函数中花费的总时间(不包括调用子函数的时间), percall 是的商 tottime 除以 ncallscumtime 是在此函数和所有子函数中花费的总时间(即,从调用到退出), percall 是的商 cumtime 除以基元调用,以及 filename:lineno(function) 提供了每个函数的相应数据。这里的经验法则是:清单中的函数越高,它的成本就越高。因此,对于优化来说,它更有趣。

像往常一样, prun? 提供有关如何使用探查器和了解输出的详细信息。

也可以将简档数据写入对象以允许更近距离的检查:

sage: %prun -r A*A
sage: stats = _
sage: stats?

注:输入 stats = prun -r A\*A 显示语法错误消息,因为prun是IPythonShell命令,而不是常规函数。

要获得良好的性能分析数据的图形表示,您可以使用HotShotProfiler,这是一个名为 hotshot2cachetree 以及该计划 kcachegrind (仅限Unix)。HotShotProfiler也是同样的示例:

sage: k,a = GF(2**8, 'a').objgen()
sage: A = Matrix(k,10,10,[k.random_element() for _ in range(10*10)])
sage: import hotshot
sage: filename = "pythongrind.prof"
sage: prof = hotshot.Profile(filename, lineevents=1)
sage: prof.run("A*A")
<hotshot.Profile instance at 0x414c11ec>
sage: prof.close()

这将生成一个文件 pythongrind.prof 在当前工作目录中。现在可以将其转换为cachegrind格式以实现可视化。

在系统Shell程序上,键入

$ hotshot2calltree -o cachegrind.out.42 pythongrind.prof

输出文件 cachegrind.out.42 现在可以使用以下命令进行检查 kcachegrind 。请注意,命名约定 cachegrind.out.XX 需要被服从。