Annotation of OpenXM/src/ox_ntl/ntl.rr, Revision 1.3
1.3 ! iwane 1: /* $OpenXM: OpenXM/src/ox_ntl/ntl.rr,v 1.2 2003/11/15 09:06:20 iwane Exp $ */
1.1 iwane 2:
3: module ntl;
4: localf factor$
1.2 iwane 5: localf lll$
1.1 iwane 6: localf ex_data$
7: localf ex_data_tmp$
1.2 iwane 8: localf mat2list$
9: localf list2mat$
1.1 iwane 10:
11:
12: /* static variables */
13:
14: /* extern variables */
15:
16: #if 1
1.3 ! iwane 17: localf test_factor$
! 18: localf test_lll$
1.1 iwane 19:
1.3 ! iwane 20: /*&usage begin: ntl.test_factor(PID, POLY)
1.1 iwane 21: compare on ox_NTL and on asir.
22:
23: example:
24: [1028] F=ntl.ex_data(4);
25: x^16-136*x^14+6476*x^12-141912*x^10+1513334*x^8-7453176*x^6+13950764*x^4-5596840*x^2+46225
26: [1029] F = F * subst(F, x, x + 1)$
27: [1030] ntl.factor(PID, F);
28: [[1,1],[x^16+16*x^15-16*x^14-1344*x^13-4080*x^12+32576*x^11+157376*x^10-255232*x^9-2062624*x^8-249088*x^7+10702080*x^6+9126912*x^5-18643712*x^4-24167424*x^3+2712576*x^2+10653696*x+2324736,1],[x^16-136*x^14+6476*x^12-141912*x^10+1513334*x^8-7453176*x^6+13950764*x^4-5596840*x^2+46225,1]]
29: [1031] ntl.test(PID, F);
30: [CPU,0.121539,0.001354,GC,0.0222,0]
31:
32: end: */
1.3 ! iwane 33: def test_factor(PID, F)
1.1 iwane 34: {
35: T0 = time();
36: fctr(F);
37: T1 = time();
38: ntl.factor(PID, F);
39: T2 = time();
40:
41: return (["CPU", T1[0]-T0[0],T2[0]-T1[0],"GC",T1[1]-T0[1],T2[1]-T1[1]]);
42: }
43:
1.3 ! iwane 44: def test_lll(PID, F)
! 45: {
! 46: T0 = time();
! 47: pari(lllint, F);
! 48: T1 = time();
! 49: ntl.lll(PID, F);
! 50: T2 = time();
! 51:
! 52: return (["CPU", T1[0]-T0[0],T2[0]-T1[0],"GC",T1[1]-T0[1],T2[1]-T1[1]]);
! 53: }
! 54:
1.1 iwane 55: #endif
56:
57: /*&usage begin: ntl.factor(PID, POLY)
58:
59: Factorize polynomial {POLY} over the rationals.
60:
61: {RETURN}
62: list
63: {POLY}
64: univariate polynomial with rational coefficients
65:
66: description:
67: Factorizes polynomial {POLY} over the rationals.
68: The result is represented by a list, whose elements are a pair represented as
69:
70: [[num,1],[factor1,multiplicity1],[factor2,multiplicity2],...].
71:
72: Products of all factor^multiplicity and num is equal to {POLY}.
73:
74: The number {num} is determined so that ({POLY}/{num}) is an integral polynomial
75: and its content (GCD of all coefficients) is 1.
76:
77: author: H.Iwane <iwane@math.sci.kobe-u.ac.jp>
78:
79: example:
80: [1282] F=(x^5-1)*(x^3+1)^2;
81: x^11+2*x^8-x^6+x^5-2*x^3-1
82: [1283] ntl.factor(PID, F);
83: [[1,1],[x^4+x^3+x^2+x+1,1],[x-1,1],[x+1,2],[x^2-x+1,2]]
84: [1284] fctr(F);
85: [[1,1],[x-1,1],[x^4+x^3+x^2+x+1,1],[x+1,2],[x^2-x+1,2]]
86:
87: ref: fctr
88: end:
89: */
90: def factor(PID, POLY)
91: {
92: local F, C, LIST, STR, RET, RET_NTL, VAR, I;
93: local TYPE;
94:
95: /* 入力チェック */
96: TYPE = type(POLY);
97: if (TYPE == 0 || TYPE == 1) {
98: return ([[POLY,1]]);
99: } else if (TYPE != 2) {
100: error("ntl.factor: invalid argument");
101: }
102:
103:
104: LIST = vars(POLY);
105: if (length(LIST) >= 2) { /* 一変数多項式のみ */
106: error("ntl.factor: invalid argument");
107: }
108:
109:
110: /* NTL で 有理係数多項式は不可 */
111: F = ptozp(POLY);
112:
113: C = sdiv(POLY, F);
114:
115: ox_cmo_rpc(PID, "fctr", F);
116:
117: RET_NTL = ox_pop_cmo(PID);
118:
119: /* ERROR Check */
120: if (type(RET_NTL) != 4 || length(RET_NTL) < 2) {
1.2 iwane 121: return (RET_NTL);
1.1 iwane 122: }
123:
124: RET = cons([RET_NTL[0] * C, 1], RET_NTL[1]);
125:
126: return (RET);
127: }
128:
129: /*&usage begin: ex_data_tmp(F, N)
130: Generate sample irreducible polynomial
131:
132: example:
133: [1032] F=t^3-p;
134: -p+t^3
135: [1033] ntl.ex_data_tmp(F,3);
136: x^27-90*x^24+1089*x^21-62130*x^18+105507*x^15-16537410*x^12-30081453*x^9-1886601330*x^6+73062900*x^3-6859000
137: [1034] ntl.ex_data_tmp(F,4);
138: -x^81+459*x^78-76896*x^75+7538094*x^72-347721147*x^69+3240161703*x^66+1032617170332*x^63-37499673798288*x^60+784360767442050*x^57+150576308695750650*x^54+771023617441694964*x^51+67248913649472410184*x^48+13913995714637027898294*x^45+270221527865987051874714*x^42-8828542741395296724347412*x^39-154971101776040822743879716*x^36+4343529580943017469231383983*x^33+4648027555241630173815780123*x^30-1072436585643253024332438894564*x^27+16237394255218510503554781142602*x^24-134104542851048701593527668875195*x^21+727430949790032393675174790142991*x^18-2727255031466780569416130788693624*x^15+7102683996190585423335589883738868*x^12-12524463688445776069953452180105904*x^9+14119870779369458313271232460210576*x^6-8591747198480108022372636571451136*x^3+2346749360904699138972190279765184
139:
140: ref: ex_data
141:
142: end: */
143: def ex_data_tmp(F, N)
144: {
145: FP = subst(F, p, prime(0));
146: GP = subst(FP, t, x);
147:
148: for (I = 1; I < N; I++) {
149: FP = subst(F, p, prime(I));
150: GP = res(t, subst(GP, x, x-t), FP);
151: }
152:
153: return (GP);
154: }
155:
156: /*&usage begin: ntl.ex_data(N)
157: Generate sample irreducible polynomial
158:
159: example:
160: [1028] ntl.ex_data(3);
161: x^8-40*x^6+352*x^4-960*x^2+576
162: [1029] ntl.ex_data(4);
163: x^16-136*x^14+6476*x^12-141912*x^10+1513334*x^8-7453176*x^6+13950764*x^4-5596840*x^2+46225
164:
165: ref:
166: ex_data_tmp
167:
168: end: */
169: def ex_data(N)
170: {
171: if (type(N) != 1 && type(N) != 10) {
172: print("invalid argument");
173: return (0);
174: }
175:
176: return (ex_data_tmp(t^2-p, N));
1.2 iwane 177: }
178:
179:
180: def mat2list(M)
181: {
182: A = size(M);
183:
184: ROW=A[0];
185: COL=A[1];
186:
187: for (I = 0; I < ROW; I++) {
188: for (J = 0; J < COL; J++) {
189: A = append(A, [M[I][J]]);
190: }
191: }
192:
193: return (A);
194: }
195:
196: def list2mat(L)
197: {
198: if (type(L) != 4) {
199: return ("Invalid Argument");
200: }
201:
202: ROW = L[0];
203: if (type(ROW) == 10)
204: ROW = int32ton(ROW);
205: COL = L[1];
206: if (type(COL) == 10)
207: COL = int32ton(COL);
208:
1.3 ! iwane 209: A = newmat(ROW, COL);
1.2 iwane 210:
211: C = 2;
212: for (I = 0; I < ROW; I++) {
213: for (J = 0; J < COL; J++) {
214: A[I][J] = L[C];
215: C++;
216: }
217: }
218:
219:
220: return (A);
221: }
222:
223:
224:
225: /*&usage begin: ntl.lll(PID, MAT)
226: the basics of LLL reducation.
227:
228: {M}
229: Matrix which element is Integer.
230:
231: example:
1.3 ! iwane 232: [1046] def trans(M) {
! 233: RET = newmat(size(M)[1], size(M)[0]);
! 234: for (I = 0; I < size(M)[0]; I++)
! 235: for (J = 0; J < size(M)[1]; J++)
! 236: RET[J][I] = M[I][J];
! 237: return (RET);
! 238: }
! 239: [1047] def lllpari(A) {
! 240: return (trans(trans(A) * pari(lllint, trans(A))));
! 241: }
! 242: [1048] M=newmat(3, 3, [[10,0, 10], [-7,3, 30], [7, 3, 20]]);
! 243: [ 10 0 10 ]
! 244: [ -7 3 30 ]
! 245: [ 7 3 20 ]
! 246: [1049] ntl.lll(PID, M);
! 247: [ 2 -6 0 ]
! 248: [ -1 -3 10 ]
! 249: [ 11 3 0 ]
! 250: [1050] lllpari(M);
! 251: [ 2 -6 0 ]
! 252: [ -1 -3 10 ]
! 253: [ 11 3 0 ]
1.2 iwane 254:
255:
256: ref:
257: pari(lll)
258:
259: end: */
260: def lll(PID, M)
261: {
262: /* parameter check */
263: TYPE = type(M);
264: if (TYPE != 6) { /* matrix */
1.3 ! iwane 265: MES = "ntl.lll: invalid argument: " + TYPE;
! 266: error(MES);
1.2 iwane 267: }
268:
269: A = mat2list(M);
270:
271: ox_cmo_rpc(PID, "lll", A);
272:
273: RET_NTL = ox_pop_cmo(PID);
274:
275: /* return value check */
276: if (type(RET_NTL) != 4) { /* list */
277: error("ntl.lll: error");
278: }
279:
280: R = list2mat(RET_NTL);
281:
282: return (R);
1.1 iwane 283: }
284:
285:
286: endmodule;
287:
288:
289: end$
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to ntl.rr CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM / src / ox_ntl