Return to taka_series.rr CVS log

Up to [local] / OpenXM / src / asir-contrib / testing

Annotation of OpenXM/src/asir-contrib/testing/taka_series.rr, Revision 1.1

1.1     ! takayama    1: /* $OpenXM$ */
        !             2:
        !             3: module tkseries;
        !             4:
        !             5: localf quote_input_form0$
        !             6: localf diff1$
        !             7: localf diffn$
        !             8: localf expand1$
        !             9: localf myeval$
        !            10:
        !            11: /* Todo */
        !            12: /* qt.dtoq, qtod;  quote<-->dp   */
        !            13: /* qt.ltoq, qtol;  quote<-->list */
        !            14:
        !            15: #define O_N 1
        !            16: #define O_P 2
        !            17: #define O_R 3
        !            18: #define O_LIST 4
        !            19: #define O_VECT 5
        !            20: #define O_MAT 6
        !            21: #define O_STR 7
        !            22: #define O_COMP 8
        !            23: #define O_DP 9
        !            24: #define O_USINT 10
        !            25: #define O_ERR 11
        !            26: #define O_GF2MAT 12
        !            27: #define O_MATHCAP 13
        !            28: #define O_F 14
        !            29: #define O_GFMMAT 15
        !            30: #define O_BYTEARRAY 16
        !            31: #define O_QUOTE 17
        !            32: #define O_OPTLIST 18
        !            33: #define O_SYMBOL 19
        !            34: #define O_RANGE 20
        !            35: #define O_TB 21
        !            36: #define O_DPV 22
        !            37: #define O_QUOTEARG 23
        !            38: #define O_VOID -1
        !            39:
        !            40: def quote_input_form0(Q) {
        !            41:   Q = quotetolist(Q);
        !            42:   return quote_input_form_quote_list(Q);
        !            43: }
        !            44: def diffn(F,X,N) {
        !            45:   if ((type(F) == O_LIST) || (type(F) == O_VECT) || (type(F) == O_MAT)) {
        !            46:     return map(diffn,F,X,N); /* It will not work for O_MAT */
        !            47:   }
        !            48:   /* We assume that F is O_QUOTE, O_P or O_R */
        !            49:   if (N == 0) return F;
        !            50:    for (I=0; I<N; I++) {
        !            51:      F = diff1(F,X);
        !            52:   }
        !            53:   return F;
        !            54: }
        !            55:
        !            56: def diff1(F,X) {
        !            57:   if (type(F) == O_QUOTE) {
        !            58:      /* It has not yet been implemented. */
        !            59:      R = "quote(calc_diff("+quote_input_form0(F)+","+rtostr(X)+"))";
        !            60:      return eval_str(R);
        !            61:   }
        !            62:   /* We assume O_P or O_R */
        !            63:   F=diff(F,x);
        !            64:   return F;
        !            65: }
        !            66:
        !            67: def myeval(C) {
        !            68:   C = objtoquote(C);
        !            69:   C = tr.simp_sin(C);  /* It is defined in test1-tr.rr */
        !            70:
        !            71:   /* Simplifier for quote should be called instead. */
        !            72:   if (qt.is_rational(C)) C = eval_quote(C);
        !            73:   return objtoquote(C);
        !            74: }
        !            75:
        !            76: /* Series expansion of F at X=C to the degree N */
        !            77: def expand1(F,X,C,N) {
        !            78:   R = 0;
        !            79:   Fn = F;
        !            80:   Cn = base_replace(Fn,[[X,C]]);
        !            81:   R = tr.simp_zero(R+myeval(Cn));
        !            82:   for (I=1; I<=N; I++) {
        !            83:      Fn = diff1(Fn,X);
        !            84:      Cn = base_replace(Fn,[[X,C]]);
        !            85:      Cn = myeval(Cn)/objtoquote(fac(I));
        !            86:      Cn = qt.cancel_number(Cn); /* simplify, e.g.,  2/4 --> 1/4 */
        !            87:      if (qt.is_minus(Cn)) Cn = qt.add_paren0(Cn); /* -1/2 --> (-1/2) */
        !            88:      Po = objtoquote(X-C)^objtoquote(I);
        !            89:      R = R + tr.simp_zero(Cn*Po);
        !            90:      R = tr.simp_zero(R);
        !            91:   }
        !            92:   Rem = eval_str("quote(oO("+rtostr(X)+","+rtostr(N+1)+"))");
        !            93:   return R+Rem;
        !            94: }
        !            95:
        !            96:
        !            97: endmodule;
        !            98:
        !            99: /* Screenshot
        !           100: [1798] load("tr.rr"); load("taka_series.rr");
        !           101: 1
        !           102:
        !           103: [1874] tkseries.expand1(cos(x),x,0,7);
        !           104: quote(1+-1/2*x^2+1/24*x^4+-1/720*x^6+oO(x, 8))
        !           105:
        !           106: By tr.simp_unary(R),  +- changes to - .
        !           107:
        !           108: [1875] print_tex_form(@@);
        !           109:      1 -  \frac{ 1} { 2}   {x}^{ 2} +  \frac{ 1} { 24}   {x}^{ 4}  -  \frac{ 1} { 720}   {x}^{ 6} + {oO}( {x},  8)
        !           110:
        !           111: [1876] tkseries.expand1(sin(x),x,0,7);
        !           112: quote(1/1*x^1+-1/6*x^3+1/120*x^5+-1/5040*x^7+oO(x, 8))
        !           113:
        !           114: [1877] print_tex_form(@@);
        !           115:       \frac{ 1} { 1}   {x}^{ 1}  -  \frac{ 1} { 6}   {x}^{ 3} +  \frac{ 1} { 120}   {x}^{ 5}  -  \frac{ 1} { 5040}   {x}^{ 7} + {oO}( {x},  8)
        !           116:
        !           117: */
        !           118:
        !           119: end$