找回密码
 欢迎注册
查看: 4242|回复: 21

[讨论] PARI/GP的一些笔记

[复制链接]
发表于 2023-11-22 16:15:04 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有账号?欢迎注册

×
在我没写完前,所有回复帖子一律删除
1、通用设置
1)提高默认堆栈容量
执行目录下有个gprc.txt
加入parisizemax = 100000000
可明显减少对大数进行复杂运算时候报堆栈错误的可能

2) 默认实数的精度
交互环境下输入\p 1000 可提高默认的实数精度,
也可以在gprc.txt里增加 realprecision = 1000
实际的精度要比输入值略大,以符合机器字的整数倍数
推荐gprc.txt默认精度不要太大,否则一次显示不全,也没必要
需要时候临时用\p 10000这种调节

3)启用计时器
交互环境下,输入#回车即可启用计时器,每次执行命令结果都会附加执行时间
再输入一次#即可关闭
##显示最后一次的执行时间

4)交互环境下,执行结果会存入%开始后跟数字的变量里
GP>x=1
%1=1
GP>y=2
%2=2
这些变量可以在输入中直接使用,
特别的单独的%代表最后一次的值

5)在线帮助
输入?会显示帮助列表,? n显示对应主题帮助
? function 会显示function函数的帮助
?? function 是更详细的function帮助
??? function 是所有和function有关联的函数清单
?\ 是显示快捷键的列表,比如#,##等
?. 是显示相关类型的成员函数


2、常见的常量
Pi 圆周率
oo 正无穷
-oo 负无穷
I 虚数单位i
Euler 欧拉常数
Catalan 卡塔兰常数

3、类型
(1)  t_INT    : Z
  无限精度整数,0x开头可以定有16进制数,0b开头定义二进制数
  例如:12345678, -234567, 0x1E2365, 0xe3, 0b10100111
  当输入整数时候,不带小数点,否则类型是实数
  GP > type(1.)
  %1 = "t_REAL"
  GP > type(1)
  %1 = "t_INT"
  
(2)  t_REAL   : R
  实数分三部分
  符号s:+1, 0, -1
  底数部分m:一个多精度整数,精度由系统变量realprecision控制
  指数部分e:根据机器字长决定的有符号整数,范围[-2^B, 2^B], 字长32位时候B=31,64时候B=63
  表示为实际值 s*m*2^e
  例如:1.0, 23.56, 3.14159, 2.3E-18, 23e29

(3)  t_INTMOD : Z/nZ
  Mod(n, m)定义n模m,其中n, m必须是整数,可以为负数,m不为0
  Mod(n, -m)=Mod(n, m)
  Mod(-n, m)=Mod(m-n, m)
  通常的+-*均可用于Mod(n, m)
  但是/,如果m不是素数可能没有结果而保存,同样倒数也可能报错

  gp > x=Mod(3,5)
  %1 = Mod(3, 5)
  gp > y=Mod(2,5)
  %2 = Mod(2, 5)
  gp > x*y
  %3 = Mod(1, 5)
  gp > x/y
  %4 = Mod(4, 5)
  gp > x+y
  %5 = Mod(0, 5)

  gp > x=Mod(4,8)
  %6 = Mod(4, 8)
  gp > y=Mod(5,8)
  %7 = Mod(5, 8)
  gp > x*y
  %8 = Mod(4, 8)
  gp > x+y
  %9 = Mod(1, 8)
  gp > x/y
  %10 = Mod(4, 8)
  gp > y/x
    ***   at top-level: y/x
    ***                  ^--
    *** _/_: impossible inverse in Fl_inv: Mod(4, 8).
  
(4)  t_FRAC   : Q
  有理数中分母不是1的形式
  分数定义为 n/m,其中n, m都是整数,且m不为0
  若m整除n,则n/m为t_INT类型
  分数运算总会化简到最简形式,这会用到GCD算法,
  通常,把分数运算转换到整数运算会提高运行速度

  gp > 2/3
  %1 = 2/3
  gp > 2^(2/3)
  %2 = 1.5874010519681994747517056392723082603914933278998530098082857618252165056242191732735442132622209570229347616813220179034976599
  gp > 27^(2/3)
  %3 = 9
  gp > 3/4+5/7
  %4 = 41/28
  gp > 1/2-1/3
  %5 = 1/6
  
(5)  t_FFELT  : Fp
  有限域类型

(6)  t_COMPLEX: C
  复数可以表示成 x + y * I,x和y可以是t_INT, t_REAL, t_FRAC, t_INTMOD, t_PADIC

  gp > x=Mod(2,65537)+Mod(1,65537)*I
  %1 = Mod(2, 65537) + Mod(1, 65537)*I
  gp > x^65537
  %2 = Mod(2, 65537) + Mod(1, 65537)*I
  gp > x^65535
  %3 = Mod(52430, 65537) + Mod(39322, 65537)*I

(7)  t_PADIC  : Qp   
  若p为素数,p_adic数定义为
  \sum_{i=k}^{oo} {a_i * p^i}
  gp > t= 2*7^(-1) + 3 + 4*7 + 2*7^2 + O(7^3)
  %1 = 2*7^-1 + 3 + 4*7 + 2*7^2 + O(7^3)
  gp > type(t)
  %2 = "t_PADIC"
  注意,系统并不检查p是不是素数,所以可以
  gp > t = 2 * (1/10 + O(10^5));
  gp > type(t)
  %4 = "t_PADIC"
  但是参与运算时候可能发生错误

(8)  t_QUAD   : Q[w]
  非平方非零整数 d, d%4=0 或者d%4=1
  w = quadgen(d, 'w)
  d%4=0时候,w=sqrt(d)/2
  d%4=1时候,w=(1+sqrt(d))/2
  二次域数表达为 x + y * w,其中x, y是t_INT, t_INTMOD,t_FRAC, t_POL,t_POLMOD类型
  令z = x + y*w
  则
  z.disc返回d
  z.mod返回z的最小多项式T
  z.pol返回t_POL类型的x+y*w
  real(z)返回x, imag(z)返回y
  gp > w=quadgen(-3,'w)
  %1 = w
  gp > z=1+2*w
  %2 = 1 + 2*w
  gp > z.disc
  %3 = -3
  gp > z.mod
  %4 = w^2 - w + 1
  gp > z.pol
  %5 = 2*w + 1
  gp > real(z)
  %6 = 1
  gp > imag(z)
  %7 = 2

(9)  t_POL    : T[X]      
  多项式,系数可以是t_INT, t_INTMOD, t_FAC, t_COMPLEX, t_QUAD

  gp > w=quadgen(-3,'w)
  %1 = w
  gp > z=1+2*w
  %2 = 1 + 2*w
  gp > p=1+x+1/2*x^2+12*x^3+Mod(1,7)*x^4+z*x^5
  %3 = (1 + 2*w)*x^5 + Mod(1, 7)*x^4 + Mod(5, 7)*x^3 + Mod(4, 7)*x^2 + Mod(1, 7)*x + Mod(1, 7)
  gp > p1=1.2+x^2
  %4 = x^2 + 1.20.....

  多项式可以用常规方法输入
  p = x^4 + 2*x^3 + 3*x^2 + 4*x + 5
  也可以用向量形式,默认变量是x,从最高次开始
  gp > P = Pol([1, 2, 3, 4, 5])
  %1 = x^4 + 2*x^3 + 3*x^2 + 4*x + 5
  也可以指定变量
  gp > P = Pol([1, 2, 3, 4, 5], 'y)
  %1 = y^4 + 2*y^3 + 3*y^2 + 4*y + 5
  也可以从最低次开始
  gp > P = Polrev([1, 2, 3, 4, 5], 'y)
  %1 = 5*y^4 + 4*y^3 + 3*y^2 + 2*y + 1

(10)  t_POLMOD : T[X]/(P)         
  MOD(x, y)中x, y都是多项式
  若z=Mod(x, y),
  z.pol返回化简后的x
  z.mod返回y
  gp > zz=Mod(x^7+1, x^2+x+1)
  %1 = Mod(x + 1, x^2 + x + 1)
  gp > zz.pol
  %2 = x + 1
  gp > zz.mod
  %3 = x^2 + x + 1
  
(11)  t_SER    : T((X))  
  类似于p-adic的形式的多项式
  gp > f = x^(-3)+2*x^(-2)+1+x^3+O(x^7)
  %1 = x^-3 + 2*x^-2 + 1 + x^3 + O(x^7)
  gp > type(f)
  %2 = "t_SER"
  非有理函数会被转化成这种形式
  gp > sin(x)
  %3 = x - 1/6*x^3 + 1/120*x^5 - 1/5040*x^7 + 1/362880*x^9 - 1/39916800*x^11 + 1/6227020800*x^13 - 1/1307674368000*x^15 + O(x^17)
  gp > type(sin(x))
  %4 = "t_SER"

(12)  t_RFRAC  : T(X)
  有理函数中分母不是1的形式
  gp > type((1+x^2)/x)
  %1 = "t_RFRAC"
  gp > type((1+x^2))
  %2 = "t_POL"

(13)  t_QFB    : 二元二次型
通过 Qfb(a, b, c) 定义的 p = a*x^2 + b*x*y + c*y^2 二次型

(14)  t_VEC    : T^n(行向量)      
(15)  t_COL    : T^n(列向量)
t_VEC和t_COL都是向量,区别是VEC是横向的,COL是纵向的
输出时候,COL有~后缀
gp > v1=Vec([1,2,3])
%1 = [1, 2, 3]                                                            ^-
gp > v2=Vec([3,2,1])
%2 = [3, 2, 1]
gp > c1=Col([1,2,3])
%3 = [1, 2, 3]~
gp > c2=Col([3,2,1])
%4 = [3, 2, 1]~
gp > m=matrix(3,3)
%5 =
[0 0 0]

[0 0 0]

[0 0 0]

gp > m[1,]=v1
%6 = [1, 2, 3]
gp > m[2,]=v2
%7 = [3, 2, 1]
gp > m[,1]=c1
%8 = [1, 2, 3]~
gp > m[,2]=c2
%9 = [3, 2, 1]~
gp > m
%10 =
[1 3 3]

[2 2 1]

[3 1 0]
两者都可以通过1开始的下标访问和赋值
v1[1], v2[3], c1[1], c2[2]等
也可以用
v[1..3]这种范围语法
v[^2]表示第二项移除后的结果
-1表示最后一个元素
因此有
V[^-1]移除最后一个元素
V[-2..-1]表示最后两个元素
也可用vector函数
vector(n, x, seq)
vector(5, i, i+1)=[2,3,4,5,6]
Col(vector(5, i, i^2)) = [1, 4, 9, 16, 25]~
或者用
[seq | x <- v, seq]
例如
[x|x<-[1..5], isprime(x)] = [2, 3, 5]
[ p | p <- primes(10), isprime(p+2) ] = [3, 5, 11, 17, 29]
[ x+y | x <- [1..3]; y <- [1..2] ] = [2, 3, 3, 4, 4, 5]
[ [x,y] | x <- [1..4], isprime(x); y <- [2..5], y % 3 == 1 ] = [[2, 4], [3, 4]]

(16)  t_MAT    : M_{m*n}(T)(矩阵)  
矩阵赋值
matrix(m,{n = m},{X},{Y},{expr = 0})
gp > m = matrix(2,2)
%1 =
[0 0]

[0 0]
gp > matrix(5, 5, x, y, x+y)
%2 =
[2 3 4 5  6]

[3 4 5 6  7]

[4 5 6 7  8]

[5 6 7 8  9]

[6 7 8 9 10]
也可以
gp > m = Mat([1,2;3,4;5,6])
%1 =
[1 2]

[3 4]

[5 6]  

m[1,], m[2,]...表示矩阵的行向量
m[,1], m[,2]...表示矩阵的列向量
m[1..2,]表示前两行
m[^1,]表示移除第一行后其他行
m[1,2..3]表示第1行,2-3列
m[2..3,1]表示第1列,2-2行
m[^1, 2..3]表示移除第1行后剩余行的2-3列

(17)  t_LIST   : 列表              
t_List类似t_VEC, t_COL,可以通过List()从他们转成t_LIST  
List() = List([])
List(Vec([1,2,3])) = List([1,2,3])
List(Col([1,2,3])) = List([1,2,3])
List([ x+y | x <- [1..3]; y <- [1..2] ]) = List([2, 3, 3, 4, 4, 5])
其他例子不再赘述
t_LIST可以通过listput追加元素
gp > l = List([1, 2, 3])
%1 = List([1, 2, 3])
gp > listput(~l, 10)
%2 = 10
gp > listput(~l, 20)
%3 = 20
gp > l
%4 = List([1, 2, 3, 10, 20])
这里的~l表示对l的引用,不需要复制到函数

(18)  t_STR    : 字符串  
  ""括起来的字符串,"用\"表示
  例如 "123"         

(19)  t_VECSMALL: 小整数向量  
  以当前机器字长限制的小整数的向量,32位下绝对值小于w^31,64位下绝对值小于2^63
  gp > x = Vec([1,2,3,4,5])
  %1 = [1, 2, 3, 4, 5]
  gp > y = Vecsmall([1,2,3,4,5])
  %2 = Vecsmall([1, 2, 3, 4, 5])
  gp > type(x)
  %3 = "t_VEC"
  gp > type(y)
  %4 = "t_VECSMALL"
  gp > x[1]=2^65
  %5 = 36893488147419103232
  gp > x
  %6 = [36893488147419103232, 2, 3, 4, 5]
  gp > y[1]=2^63
  ***   at top-level: y[1]=2^63
  ***                 ^---------
  ***   incorrect type in t_VECSMALL assignment (t_INT).
  
(20)  t_CLOSURE: 函数
  所有函数

(21)  t_ERROR  : 错误文本     
每当发生错误时都会创建这种类型的对象:它包含有关错误和错误上下文的一些信息。通常,会立即打印适当的错误,计算中止,GP进入“中断循环”
gp > 1/0; 1 + 1
***   at top-level: 1/0;1+1
***                  ^------
*** _/_: division by a noninvertible object
***   Break loop: type 'break' to go back to the GP prompt

可以用iferr捕获错误,并打印错误信息
gp > i = Mod(6,12); iferr(1/i, E, print(E)); 1 + 1
error("impossible inverse modulo: Mod(6, 12).")
%1 = 2
gp > i = Mod(6,12); iferr(1/i, E, print(component(E,1)));
Mod(6, 12)

(22)  t_INFINITY: -oo 和 +oo
  正负无穷      

点评

费了不少功夫重新补充了一大堆内容~~~主要是类型  发表于 2023-11-24 15:22
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2023-11-22 16:44:59 | 显示全部楼层
4、函数定义和常见函数
1)定义

gp > f(x)=x^2+1
%1 = (x)->x^2+1
gp > f(3)
%2 = 10

gp > f(x,y)=x+y
%3 = (x,y)->x+y
gp > f(i,5)
%4 = i + 5

或者以下方式也可以
f=x->x^2+1
f=(x)->x^2+1
f=(x,y)->x+y

也可以定义多行函数,每行用;结束即可,函数值是最后一行的结果

g(x) = y = x^2; z = y % 9;
g(7)
  4

如果在某些需要seq参数的函数中,这些行需要用{}括起来


点评

nyy
弄个PDF或者  发表于 2023-11-22 19:39
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2023-11-22 17:05:22 | 显示全部楼层
2) if
if(a, seq1)
如果 a != 0 执行seq1
if(a, seq1, seq2)
如果 a != 0 执行seq1,否则执行seq2
if(a1, seq1,
   a2, seq2.
   .....
    defaultseq)
这个执行多个判定条件,类似连续的if判断
所有条件都不满足的话,执行最后的默认seq
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2023-11-22 17:47:04 | 显示全部楼层
3) for
for(x=a, b, seq)
a, b是实数,从a 到 b 执行 seq, 步长是1
如果 a > b 则不执行
如果b 是 +∞(b = oo ),则循环是无限循环

GP > for(x=2.5, 3, print(x))
2.5

GP > for(x=2.5, 4, print(x))
2.5
3.5

可以用break跳出循环
GP> for(i=1,100,{print(i); if(i>=10,break())})
1
2
3
4
5
6
7
8
9
10
可以用next跳过当前循环进入下一个
gp > for(i=1,10,if(isprime(i), next, print(i)))
1
4
6
8
9
10

其他for函数:
forcomposite(n = a, b, seq)
对[a, b] 区间的合数执行 seq
gp > forcomposite(a=2, 10,print(a))
4
6
8
9
10

forprime(n = a, b, seq)
对[a, b] 区间的素数执行 seq
gp > forprime(n=2,10, print(n))
2
3
5
7

forprimestep(p = a, b, q, seq)
这里a, b都是unsigned long类型(32或者64位), p是素数, a <= p <= b
若q是t_INT,p满足p=a mod q
若q是t_INTMOD, q = Mod(c, N), 则Mod(p, N)=Mod(c, N)
gp > forprimestep(p = 4, 30, 5, print(p))
   19
   29
gp > forprimestep(p = 4, 30, Mod(1,5), print(p))
   11

fordiv(n, d, seq)
对n的因子d执行 seq
gp > fordiv(100, d, print(d))
1
2
4
5
10
20
25
50
100

fordivfactored(n, d, seq)
对n进行分解,然后对n的所有因子和因子的分解执行seq

gp > fordivfactored(24, d, print(d))
[1, matrix(0,2)]
[2, Mat([2, 1])]
[3, Mat([3, 1])]
[4, Mat([2, 2])]
[6, [2, 1; 3, 1]]
[8, Mat([2, 3])]
[12, [2, 2; 3, 1]]
[24, [2, 3; 3, 1]]
gp > fordivfactored(24, d, print(eulerphi(d)))
1
1
2
2
2
4
4
8

forfactored(n = a, b, seq)
从a到b,对n和n的分解执行seq
gp > forfactored(n=2,10, print(n))
[2, Mat([2, 1])]
[3, Mat([3, 1])]
[4, Mat([2, 2])]
[5, Mat([5, 1])]
[6, [2, 1; 3, 1]]
[7, Mat([7, 1])]
[8, Mat([2, 3])]
[9, Mat([3, 2])]
[10, [2, 1; 5, 1]]


forsquarefree(N=a, b, seq)
无平方因子的N和N的分解[N, factor(N)]
gp > s=0;forsquarefree(n=1,100, s+=n[1]);s
%1 = 2967
gp > forsquarefree(n=1,10, print(n));
[1, matrix(0,2)]
[2, Mat([2, 1])]
[3, Mat([3, 1])]
[5, Mat([5, 1])]
[6, [2, 1; 3, 1]]
[7, Mat([7, 1])]
[10, [2, 1; 5, 1]]

foreach(v, x, seq)
对v中每个元素x执行seq, v的类型是t_LIST, t_VEC,t_COL, t_MAT
如果类型是t_MAT,x是矩阵的行
gp > foreach(factor(24), d, print(d))
[2, 3]~
[3, 1]~


forpart(X = k,seq,{a = k1},{n = k2})
X是满足x_1+x_2...+x_n = k的vec,如果可选项a=k1存在,则x_i最大是k1,如果可选项n=k2存在,则最多k2项
gp > forpart(X=4, print(X))
   Vecsmall([4])
   Vecsmall([1, 3])
   Vecsmall([2, 2])
   Vecsmall([1, 1, 2])
   Vecsmall([1, 1, 1, 1])
gp > forpart(v=5,print(Vec(v)), 4, 3)
   [1, 4]
   [2, 3]
   [1, 1, 3]
   [1, 2, 2]
有2-4项,小于5,填充0
gp > forpart(v=5,print(Vec(v)),[0,5],[2,4])
   [0, 0, 1, 4]
   [0, 0, 2, 3]
   [0, 1, 1, 3]
   [0, 1, 2, 2]
   [1, 1, 1, 2]

forperm(a, p, seq)
若a是t_INT,对1..a的全排列执行seq
若a是t_VECS,对a的全排列执行seq
gp > forperm(3, p, print(p))
   Vecsmall([1, 2, 3])
   Vecsmall([1, 3, 2])
   Vecsmall([2, 1, 3])
   Vecsmall([2, 3, 1])
   Vecsmall([3, 1, 2])
   Vecsmall([3, 2, 1])
gp > forperm([2,1,1,3], p, print(p))
   Vecsmall([2, 1, 1, 3])
   Vecsmall([2, 1, 3, 1])
   Vecsmall([2, 3, 1, 1])
   Vecsmall([3, 1, 1, 2])
   Vecsmall([3, 1, 2, 1])
   Vecsmall([3, 2, 1, 1])

forstep(X = a,b,s,seq)
a<=X<=b
若s是整数,则X=a, a+s, a+2s,a+3s......
若s是Mod(c, N),则Mod(X, N) = Mod(c, N)
若s是vec, 比如s = [1, 2]
则X = a, a+1, a+1+2, a+1+2+1, a+1+2+1+2...
gp > forstep(x=5, 10, 2, print(x))
   5
   7
   9
gp > forstep(x=5, 10, Mod(1,3), print(x))
   7
   10
gp > forstep(x=5, 10, [1,2], print(x))
   5
   6
   8
   9
gp >  forstep(x=5, 20, [1,2, 3], print(x))
5
6
8
11
12
14
17
18
20

forsubgroup(H = G,{bound},seq)
群G的阿贝尔子群,执行seq

forsubset(nk, s, seq)
若nk是整数,[1,2,3.....,nk]的所有子集执行seq
若nk=[n, k], [1,2,3.....,nk]的所有k个元素的子集执行seq
gp > forsubset([5,3], s, print(s))
   Vecsmall([1, 2, 3])
   Vecsmall([1, 2, 4])
   Vecsmall([1, 2, 5])
   Vecsmall([1, 3, 4])
   Vecsmall([1, 3, 5])
   Vecsmall([1, 4, 5])
   Vecsmall([2, 3, 4])
   Vecsmall([2, 3, 5])
   Vecsmall([2, 4, 5])
   Vecsmall([3, 4, 5])
gp > forsubset(3, s, print(s))
   Vecsmall([])
   Vecsmall([1])
   Vecsmall([2])
   Vecsmall([3])
   Vecsmall([1, 2])
   Vecsmall([1, 3])
   Vecsmall([2, 3])
   Vecsmall([1, 2, 3])

forvec(X = v,seq,{flag = 0})
v是n个[a_i, b_i]组成的,
函数生成n个元素的vec,每个元素v_i,都满足a_i<=v_i<=b_i
如果flag=1,vec中元素满足v_i<=v_j, i<j
如果flag=2,vec中元素满足v_i<v_j, i<j

gp > forvec (X=[[0,1],[-1,1]], print(X));
   [0, -1]
   [0, 0]
   [0, 1]
   [1, -1]
   [1, 0]
   [1, 1]
gp > forvec (X=[[0,1],[-1,1]], print(X), 1);
   [0, 0]
   [0, 1]
   [1, 1]
   ? forvec (X=[[0,1],[-1,1]], print(X), 2)
   [0, 1]

gp >  forvec (X=[[0,2],[-1,1], [2, 4]], print(X));
[0, -1, 2]
[0, -1, 3]
[0, -1, 4]
[0, 0, 2]
[0, 0, 3]
[0, 0, 4]
[0, 1, 2]
[0, 1, 3]
[0, 1, 4]
[1, -1, 2]
[1, -1, 3]
[1, -1, 4]
[1, 0, 2]
[1, 0, 3]
[1, 0, 4]
[1, 1, 2]
[1, 1, 3]
[1, 1, 4]
[2, -1, 2]
[2, -1, 3]
[2, -1, 4]
[2, 0, 2]
[2, 0, 3]
[2, 0, 4]
[2, 1, 2]
[2, 1, 3]
[2, 1, 4]
gp >  forvec (X=[[0,2],[-1,1], [2, 4]], print(X), 1);
[0, 0, 2]
[0, 0, 3]
[0, 0, 4]
[0, 1, 2]
[0, 1, 3]
[0, 1, 4]
[1, 1, 2]
[1, 1, 3]
[1, 1, 4]
gp >  forvec (X=[[0,2],[-1,1], [2, 4]], print(X), 2);
[0, 1, 2]
[0, 1, 3]
[0, 1, 4]

毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2023-11-23 08:39:52 | 显示全部楼层
4)while
while(a, seq)
当a不为0时候,执行seq,直到a是0

gp > x=1;while(x<=10, {print(x);x=x+1});
1
2
3
4
5
6
7
8
9
10

注意,
如果a永远为1,将无限循环seq
如果a一开始为0,则不执行seq
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2023-11-23 08:47:02 | 显示全部楼层
5) until
until(a, seq)
执行seq直到a为1
需要注意的是,seq会至少执行一次

gp > x=1;until(x>0, {print(x);x=x+1});
1
gp > x=1;until(x>10, {print(x);x=x+1});
1
2
3
4
5
6
7
8
9
10

until 和 while 是相反的条件
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2023-11-27 14:35:37 | 显示全部楼层
6)输入输出

(1)input(): 获取命令行输入
gp > x = input()
1
%1 = 1

(2)print({str}*),print1({str}*), printp({str}*): 输出
gp > x = 1
%1 = 1
gp > print(x)
1
gp > m = [1,2;3,4]
%83 =
[1 2]

[3 4]

gp > print(m) //以raw方式显示,等价于print1
[1, 2; 3, 4]
gp > printp(m) //prettymatrix方式显示

[1 2]

[3 4]

(3) printsep(sep,{str}*)
用sep字符间隔的输出
printsep1(sep,{str}*)
和printsep类似,但是输出成raw形式

gp > printsep(":", 1,2,3,4)
1:2:3:4
gp > printsep(".", 1,2,3,4)
1.2.3.4
gp > printsep(",", 1,2,3,4)
1,2,3,4
gp > printsep(" ", 1,2,3,4)
1 2 3 4
gp > printsep1(":", [1,2;3,4], [5,6;7,8])
[1, 2; 3, 4]:[5, 6; 7, 8]

(4) printf(fmt, {x}*)
根据fmt字符串格式化输出,类似C的同名函数
不带%的原样输出
%的根据转换规则输出
即%开始,
    后跟一个或多个可选标注(#, 0, -, +, " "),
    一个可选的输出最小输出字符串长度,一个可选的小数精度,这两个用 . 分割,
    可选的转换类型字符 (l, d, i, o, u, x, X, p, e, E, f, F,g, G, s)
#表示16进制输出x格式加前缀"0x", X格式加前缀"0X",8进制o格式加前缀"o"
0表示数字类型的,规定最小输出字符串长度后,不足部分补0,
" "表示,规定最小输出字符串长度后,不足部分补空格
-表示输出会对齐在左边,默认对齐右边,- 和 0同时出现,忽略0
        即当输出长度大于实际数字宽度时候,- 会导致在输出位置先输出数字,后出现补齐的空格
+表示 + 和 - 符号总是会输出,+ 和" "同时出现,会忽略" "
m.n 宽度限定,表示最少输出m个字符,小数固定n位,m可省略
类型部分
l: 有符号机器字十进制整数
d, i:有符号十进制整数
u:无符号十进制整数
(当用u输出大的负整数时候,显示的是补码表示)
o: 无符号八进制整数
x, X:无符号十六进制表示,其中A-F字母,x格式用小写a-f表示,X格式用大写A-F表示
p: GEN类型
e, E:数字表达成小数 d * 10^n 形式,其中1<= d < 10, n是有符号整数, e格式输出成d en, E格式输出成d En
f, F:当数字绝对值小于1时候,总是输出成小数 d.dddd 的形式,
        当数字绝对值大于等于1时候,如果小数形式精度可以输出到小数点,用小数形式,
        否则输出和e, E一致(这里要依赖于realprecision的设定)
                例如realprecision = 19, 输出 123456789012345678.9可以用小数形式,但1234567890123456789就得用指数形式
g, G: 当输出数字绝对值小于1,如果绝对值小于1E-4,用小数形式输出,否则等同于e, E
        当输出数字绝对值大于等于1,按照 f, F形式处理
s: 字符串

gp > printf("floor: %d, field width 3: %3d, with sign: %+3d\n", Pi, 1, 2);
   floor: 3, field width 3:   1, with sign:  +2

gp > printf("%.5g %.5g %.5g\n",123,123/456,123456789);
   123.00 0.26974 1.2346 e8

gp >  x = 23; y=-1/x; printf("x=%+06.2f y=%+0*.*f\n", x, 6, 2, y);
   x=+23.00 y=-00.04

gp > for (i = 2, 5, printf("%05d\n", 10^i))
   00100
   01000
   10000
   100000  \\ don't truncate fields whose length is larger than the minimum width
gp > printf("%.2f  |%06.2f|", Pi,Pi)
   3.14  |  3.14|

gp > printf("%4d", [1,2,3]);
   [   1,   2,   3]
gp >  printf("%5.2f", mathilbert(3));
   [ 1.00  0.50  0.33]

   [ 0.50  0.33  0.25]

   [ 0.33  0.25  0.20]

(5) strprintf(fmt, {x}*)
和printf类似,但是会输出到一个字符串里

(6) printtex({x}*), strtex({x}*): 输出成tex形式
printtex在标准输出设备上输出,strtex输出到字符串
gp >  v = [1, 2, 3]
   %1 [1, 2, 3]
   ? strtex(v)
   %2 = "\\pmatrix{ 1&2&3\\cr}\n"
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2023-11-27 16:37:55 | 显示全部楼层
7) strjoin 和 strsplit
strjoin(v,{p = ""}):
  把v中元素用p连接起成字符串,如果p未指定,则v直接连接成字符串
gp > strjoin([1,2,3,4,5], "-")
%1 = "1-2-3-4-5"
gp > strjoin([1,2,3,4,5], "*")
%2 = "1*2*3*4*5"
gp > strjoin([1,2,3,4,5], ".")
%3 = "1.2.3.4.5"
gp > strjoin([1,2,3,4,5])
%4 = "12345"
gp > strjoin(["hello", "world"], ",")
%5 = "hello,world"


strsplit(s,{p = ""}):
和strjoin相反的操作,把字符串 s 通过 p,分割成 vec, 如果 p 未指定,则分割成单个字符序列
gp > strsplit("1-2-3-4-5", "-")
%1 = ["1", "2", "3", "4", "5"]
gp > strsplit("1-2-3-4-5")
%1 = ["1", "-", "2", "-", "3", "-", "4", "-", "5"]
如果开始,结束位置有p字符,将分割出0长度字符""
gp > strsplit("-1-2-3-4-5-", "-")
%1 = ["", "1", "2", "3", "4", "5", ""]
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2023-11-27 16:40:00 | 显示全部楼层
8) apply, select

apply(f, A)
如果A是数字或者字符串,返回f(A)
如果A是t_POL或者t_SER,返回对应类型,但是对应系数a_i系数变成f(a_i)
如果A是t_VEC, t_VECSMALL, t_LIST, t_MAT,对应位置的值x变成f(x)
gp > f=x->x^2+1
%1 = (x)->x^2+1
gp > apply(f, 2)
%2 = 5
gp > apply(f, [1,2,3])
%3 = [2, 5, 10]
gp > apply(f, y^2+2*y+1)
%4 = 2*y^2 + 5*y + 2
gp > apply(f, [1,2;3,4])
%5 =
[ 2  5]

[10 17]

select(f, A, {flag = 0})、
A是t_VEC 或者 t_LIST
如果flag=0,返回 A 中满足 f(x) 不是0 的元素组成的新 t_VEC 或者 t_LIST
如果flag=1,返回 A 中满足 f(x) 不是0 的元素在原变量里的位置索引组成的 t_VECSMALL
gp > select(isprime, vector(100, i, i))
%1 = [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97]
gp > select(isprime, vector(100, i, i^2+1))
%2 = [2, 5, 17, 37, 101, 197, 257, 401, 577, 677, 1297, 1601, 2917, 3137, 4357, 5477, 7057, 8101, 8837]
gp > select(isprime, vector(100, i, i^2+1), 1)
%3 = Vecsmall([1, 2, 4, 6, 10, 14, 16, 20, 24, 26, 36, 40, 54, 56, 66, 74, 84, 90, 94])
gp > select(isprime, List([2,3,4,5,6,7,8,9]))
%4 = List([2, 3, 5, 7])
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2023-11-27 16:56:22 | 显示全部楼层
9)fold, prod, sum

fold(f, A)
f是二元函数f(x, y),A是至少有两个元素的t_VEC
返回f(...f(f(A[1],A[2]),A[3])...,A[#A])
gp > f=(x,y)->x+y
%1 = (x,y)->x+y
gp > fold(f, [1..5])
%2 = 15

prod(X=a, b, expr, {x=1})
遍历a, b然后计算expr, 再把expr结果乘起来
如果最后的参数x默认不是1,再乘以x
gp > prod(i=1, 10, i)
%1 = 3628800
gp > prod(i=1, 10, i, 10)
%2 = 36288000

sum(X = a,b,expr,{x = 0}):
遍历a, b然后计算expr, 再把expr结果加起来
如果最后的参数x默认不是0,再加x
gp > sum(i=1, 10, i)
%1 = 55
gp > sum(i=1, 10, i, 10)
%2 = 65

点评

@nyy, 这个事你来做  发表于 2023-11-28 10:44
nyy
出个pdf  发表于 2023-11-27 20:56
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
您需要登录后才可以回帖 登录 | 欢迎注册

本版积分规则

小黑屋|手机版|数学研发网 ( 苏ICP备07505100号 )

GMT+8, 2024-12-4 05:00 , Processed in 0.036132 second(s), 17 queries .

Powered by Discuz! X3.5

© 2001-2024 Discuz! Team.

快速回复 返回顶部 返回列表