[BACK]Return to ntl.rr CVS log [TXT][DIR] Up to [local] / OpenXM / src / ox_ntl

Annotation of OpenXM/src/ox_ntl/ntl.rr, Revision 1.4

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

FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>