[BACK]Return to bfunction.sm1 CVS log [TXT][DIR] Up to [local] / OpenXM / src / kan96xx / Doc

Annotation of OpenXM/src/kan96xx/Doc/bfunction.sm1, Revision 1.2

1.1       maekawa     1: %%% oaku/kan/bfunction.sm1, 1998, 11/5
                      2:
                      3: %%% global variables for bfunction
                      4: %%% bfunction.*
                      5: /bfunction.version (2.981105) def
                      6: bfunction.version [(Version)] system_variable gt
                      7: { (This package requires the latest version of kan/sm1) message
                      8:   (Please get it from http://www.math.kobe-u.ac.jp/KAN) message
                      9:   error
                     10: } { } ifelse
                     11: /bfunction.v [(x) (y) (z)] def   %% default variables of the input polynomial
                     12: /bfunction.s (s) def             %% default variable of the output b-function
                     13: /bfunction.vh (v_) def           %% default variable for V-homogenization
                     14: /bfunction.t (t_) def            %% default variable for t in delta(t-f)
                     15: /bfunction.a  [] def             %% parameters are not available yet
                     16: /bfunction.verbose 0 def         %% no messages if 0
                     17: /bfunction.strategy 0 def        %% V-homogenization + h-homogenization if 0
                     18:                                  %% V-homogenization if 1 (not available yet)
                     19:                                  %% h-homogenization if 2 (not available yet)
                     20: /bfunction.result 0 def
                     21:
                     22: (bfunction.sm1,  11/05,1998  (C) T. Oaku.  bfunction ) message-quiet
                     23:
                     24: [(bfunction)
                     25:  [( a bfunction b)
                     26:   (array a; poly b;)
                     27:   (a :  [f] ;  string f ;)
                     28:   (a :  [f] ;  polynomial f ;)
                     29:   (a :  [f v] ; string f,v; )
                     30:   (a :  [f v] ; polynomial f, string v; )
                     31:   (b is the b-function (=Bernstein-Sato polynomial) of a polynomial f)
                     32:   (in variables v.)
                     33:   (If v is not specified, the variables are assumed to be (x,y,z). )
                     34:   (b will be a polynomial in s.  This variable can be changed by typing in)
                     35:   (  (variable) /bfunction.s set )
                     36:   (For the algorithm, see Duke Math. J. 87 (1997),115-132,)
                     37:   (   J. Pure and Applied Algebra 117&118(1997), 495--518.)
                     38:   $Example  [(x^3-y^2) (x,y)] bfunction :: $
1.2     ! takayama   39:   (  )
        !            40:   (See also bfct which implements a new algorithm to compute b-function and is faster. Aug 2002.)
1.1       maekawa    41:  ]
                     42: ] putUsages
                     43:
                     44: /bfunction {
                     45:   /arg1 set
                     46:   [/aa /typev /setarg /f /s /v /bf /bfs /vt ] pushVariables
                     47:   [(CurrentRingp) (KanGBmessage)] pushEnv  %% push current global environment.
                     48:   [
                     49:
                     50:     /aa arg1 def
                     51:     aa isArray { } { (array bfunction) message error } ifelse
                     52:     /setarg 0 def
                     53:     aa { tag } map /typev set
                     54:     typev [ StringP ] eq
                     55:     {  /f aa 0 get def
                     56:        /v bfunction.v def
                     57:        /s bfunction.s def
                     58:        /setarg 1 def
                     59:     } { } ifelse
                     60:     typev [ PolyP ] eq
                     61:     {  /f aa 0 get (string) data_conversion def
                     62:        /v bfunction.v def
                     63:        /s bfunction.s def
                     64:        /setarg 1 def
                     65:     } { } ifelse
                     66:     typev [StringP StringP] eq
                     67:     {  /f aa 0 get def
                     68:        /v [ aa 1 get to_records pop ] def
                     69:        /s bfunction.s def
                     70:        /setarg 1 def
                     71:     } { } ifelse
                     72:     typev [PolyP StringP] eq
                     73:     {  /f aa 0 get (string) data_conversion def
                     74:        /v [ aa 1 get to_records pop ] def
                     75:        /s bfunction.s def
                     76:        /setarg 1 def
                     77:     } { } ifelse
                     78:     typev [StringP ArrayP] eq
                     79:     {  /f aa 0 get def
                     80:        /v aa 1 get def
                     81:        /s bfunction.s def
                     82:        /setarg 1 def
                     83:     } { } ifelse
                     84:     typev [PolyP ArrayP] eq
                     85:     {  /f aa 0 get (string) data_conversion def
                     86:        /v aa 1 get def
                     87:        /s bfunction.s def
                     88:        /setarg 1 def
                     89:     } { } ifelse
                     90:     setarg { } { (Argument mismatch) message error } ifelse
                     91:
                     92:     [(KanGBmessage) bfunction.verbose] system_variable
                     93:
                     94:     v bfunction.t append /vt set
                     95:
                     96:     [f v fw_delta bfunction.t vt] indicial 0 get /bf set
                     97:     [bfunction.s ring_of_polynomials 0] define_ring
                     98:     bf . /bf set
                     99:     bfunction.s . /bfs set
                    100:     bf [[bfs (-1). bfs sub]] replace /bf set
                    101:     /bfunction.result bf def
                    102:     /arg1 bf def
                    103:   ] pop
                    104:   popEnv
                    105:   popVariables
                    106:   arg1
                    107: } def
                    108:
                    109: %%  Computing the indicial polynomial (the b-function) of a D-module
                    110: /indicial {
                    111:   /arg1 set  %% [equations, the variable to be restricted to 0, all variables]
                    112:   [/eqs /t /vars /allvars /newvars /x_vars /ans1 /ans2 ] pushVariables
                    113:   [(CurrentRingp)] pushEnv
                    114:   [
                    115:     arg1 0 get /eqs set
                    116:     arg1 1 get /t set
                    117:     arg1 2 get /vars set
                    118:     vars bfunction.s append /allvars set
                    119:     [bfunction.t] allvars complement /newvars set
                    120:     [bfunction.t] vars complement /x_vars set
                    121:     [eqs t vars] indicial1 /ans1 set
                    122:     [ans1 x_vars newvars] eliminate_Dx /ans2 set
                    123:     [ans2 x_vars newvars] eliminate_x /arg1 set
                    124:   ] pop
                    125:   popEnv
                    126:   popVariables
                    127:   arg1
                    128: } def
                    129:
                    130: %% (-1,0;1,0)-Groebner basis
                    131: %% [equations (t) vars] indical1 ---> psi(BFequations) (as a list of strings)
                    132: /indicial1 {
                    133:  /arg1 set
                    134:  [/bft /bfs /bfh /bf1 /ff /ans /n /i /BFallvarlist /BFDvarlist
                    135:   /BFs_weight /BFvarlist ] pushVariables
                    136:  [(CurrentRingp)] pushEnv
                    137:  [
                    138:   /ff arg1 0 get def
                    139:   /bft arg1 1 get def
                    140:   /BFvarlist arg1 2 get def
                    141:   /BFallvarlist
                    142:     [ bfunction.vh bfunction.s] BFvarlist concat bfunction.a concat
                    143:   def
                    144:   BFvarlist length /n set
                    145:   BFvarlist {xtoDx} map /BFDvarlist set
                    146:   /BFs_weight
                    147:     [ [ bfunction.vh 1 ]
                    148:       [ 0 1 n 1 sub
                    149:           { /i set BFDvarlist i get 1 }
                    150:         for
                    151:         0 1 n 1 sub
                    152:           { /i set BFvarlist i get 1 }
                    153:         for ]
                    154:     ] def
                    155:
                    156:   [ BFallvarlist listtostring ring_of_differential_operators
                    157:     BFs_weight weight_vector
                    158:   0] define_ring
                    159:
                    160:   /bfh  (h). def
                    161:   /bfs  bfunction.vh . def
                    162:   /bf1  (1). def
                    163:   ff { bft fw_homogenize . } map /ff set
                    164:   ff {[[bfh bf1]] replace} map {homogenize} map /ff set
                    165:   [ff] groebner 0 get {[[bfh bf1]] replace} map /ff set
                    166:   ff reducedBase /ans set
                    167:   ans {bft fw_psi} map /ans set
                    168:   ans {(string) data_conversion} map /arg1 set
                    169:   ] pop
                    170:   popEnv
                    171:  popVariables
                    172:  arg1
                    173: } def
                    174:
                    175: %% eliminates Dx in the ring of differential operators
                    176: /eliminate_Dx {
                    177:  /arg1 set  %% [operators x variables]
                    178:  [/bfh /bf1 /ff /ans /nx /ny /x_varlist /Dx_weight /BFvarlist
                    179:    /allvarlist /Dx_varlist /y_varlist /Dy_varlist /allvarlist /i
                    180:  ] pushVariables
                    181:  [(CurrentRingp)] pushEnv
                    182:  [
                    183:   /ff arg1 0 get def
                    184:   /x_varlist arg1 1 get def
                    185:   /BFvarlist arg1 2 get def
                    186:   x_varlist length /nx set
                    187:   BFvarlist bfunction.a concat /allvarlist set
                    188:
                    189:   x_varlist {xtoDx} map /Dx_varlist set
                    190:   x_varlist BFvarlist complement /y_varlist set
                    191:   y_varlist length /ny set
                    192:   y_varlist {xtoDx} map /Dy_varlist set
                    193:
                    194:   /Dx_weight
                    195:     [ [ 0 1 nx 1 sub
                    196:           { /i set Dx_varlist i get 1 }
                    197:         for ]
                    198:       [ 0 1 nx 1 sub
                    199:           { /i set x_varlist i get 1 }
                    200:         for
                    201:         0 1 ny 1 sub
                    202:           { /i set y_varlist i get 1 }
                    203:         for
                    204:         0 1 ny 1 sub
                    205:           { /i set Dy_varlist i get 1 }
                    206:         for
                    207:       ]
                    208:     ] def
                    209:
                    210:   [ allvarlist listtostring  ring_of_differential_operators
                    211:     Dx_weight weight_vector
                    212:   0] define_ring
                    213:
                    214:   /bfh (h). def
                    215:   /bf1 (1). def
                    216:   ff {.} map /ff set
                    217:   ff {[[bfh bf1]] replace} map {homogenize} map /ff set
                    218:   bfunction.verbose 1 eq
                    219:     {(Eliminating the derivations w.r.t. ) messagen x_varlist ::}
                    220:     { }
                    221:   ifelse
                    222:   [ff] groebner 0 get {[[bfh bf1]] replace} map /ff set
                    223:   ff reducedBase /ans set
                    224:   ans Dx_varlist eliminatev /ans set
                    225:   ans {(string) data_conversion} map /arg1 set
                    226:  ] pop
                    227:  popEnv
                    228:  popVariables
                    229:  arg1
                    230: } def
                    231:
                    232: %% eliminates x in the ring of polynomials
                    233: /eliminate_x {
                    234:  /arg1 set  %% [operators x variables]
                    235:  [/bfh /bfs /bf1 /ff /ans /nx /ny /x_varlist  /BFvarlist
                    236:    /allvarlist /y_varlist /i
                    237:  ] pushVariables
                    238:  [(CurrentRingp)] pushEnv
                    239:  [
                    240:   /ff arg1 0 get def
                    241:   /x_varlist arg1 1 get def
                    242:   /BFvarlist arg1 2 get def
                    243:   x_varlist length /nx set
                    244:   BFvarlist bfunction.a concat /allvarlist set
                    245:
                    246:   x_varlist BFvarlist complement /y_varlist set
                    247:   y_varlist length /ny set
                    248:
                    249:   /x_weight
                    250:     [ [ 0 1 nx 1 sub
                    251:           { /i set x_varlist i get 1 }
                    252:         for ]
                    253:       [ 0 1 ny 1 sub
                    254:           { /i set y_varlist i get 1 }
                    255:         for
                    256:       ]
                    257:     ] def
                    258:
                    259:   [ allvarlist listtostring  ring_of_polynomials x_weight weight_vector
                    260:   0] define_ring
                    261:
                    262:   /bfh (h). def
                    263:   /bf1 (1). def
                    264:   ff {.} map /ff set
                    265:   ff {[[bfh bf1]] replace} map {homogenize} map /ff set
                    266:   bfunction.verbose 1 eq
                    267:     {(Eliminating the variables ) messagen x_varlist ::}
                    268:     { }
                    269:   ifelse
                    270:   [ff] groebner 0 get {[[bfh bf1]] replace} map /ff set
                    271:   ff reducedBase /ans set
                    272:   ans x_varlist eliminatev /ans set
                    273:   ans {(string) data_conversion} map /arg1 set
                    274:  ] pop
                    275:  popEnv
                    276:  popVariables
                    277:  arg1
                    278: } def
                    279: %%%%%%%%%%%%%%%%%%%%%%% libraries  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                    280:
                    281: %% FW-principal part of an operator (FW-homogeneous)
                    282: %%  Op (poly) fw_symbol --->  FW-symbol(Op)  (poly)
                    283: /fw_symbol {
                    284:   [[(h). (1).]] replace bfunction.vh . coefficients 1 get 0 get
                    285: } def
                    286:
                    287: %% FW-homogenization
                    288: %% Op (string) (t) fw_homogenize ---> h(Op) (string)
                    289: /fw_homogenize {
                    290:   /arg2 set  %% bft (string)
                    291:   /arg1 set  %% an operator (string)
                    292:   [ /bft /bfDt /bfht /bfhDt /Op /degs /m /mn ] pushVariables
                    293:   [
                    294:     /Op arg1 expand def
                    295:     /bft arg2 def
                    296:     bft xtoDx /bfDt set
                    297:     bfunction.vh (^(-1)*) bft 3 cat_n /bfht set
                    298:     bfunction.vh (*) bfDt 3 cat_n /bfhDt set
                    299:     Op [[bft expand bfht expand][bfDt expand bfhDt expand]] replace
                    300:       /Op set
                    301:     Op bfunction.vh expand coefficients 0 get
                    302:       {(integer) data_conversion} map /degs set
                    303:     degs << degs length 1 sub >> get /m set
                    304:     0 m sub /mn set
                    305:     << bfunction.vh expand mn powerZ >> Op mul /Op set
                    306:     Op (string) data_conversion /arg1 set
                    307:   ] pop
                    308:   popVariables
                    309:   arg1
                    310: } def
                    311:
                    312: %% setup the ring of differential operators with the variables varlist
                    313: %% and parameters bfunction.a
                    314: %% varlist setupBFring
                    315: /setupDring {
                    316:   /arg1 set
                    317:   [ /varlist /bft /allvarlist /n /dvarlist /D_weight /i
                    318:   ] pushVariables
                    319:   [
                    320:     arg1 /varlist set
                    321:     /allvarlist
                    322:       varlist bfunction.a join
                    323:     def
                    324:     varlist length /n set
                    325:     varlist {xtoDx} map /dvarlist set
                    326:     /D_weight
                    327:     [ [ 0 1 n 1 sub
                    328:           { /i set dvarlist i get 1 }
                    329:         for ]
                    330:       [
                    331:         0 1 n 1 sub
                    332:           { /i set varlist i get 1 }
                    333:         for ]
                    334:     ] def
                    335:
                    336:     [ allvarlist listtostring ring_of_differential_operators
                    337:       D_weight weight_vector
                    338:     0] define_ring
                    339:
                    340:   ] pop
                    341:   popVariables
                    342: } def
                    343:
                    344: %% psi(P)(s)
                    345: %% Op (poly) (t) (string) fw_psi ---> psi(P) (poly)
                    346: %% Op should be FW-homogeneous.
                    347: /fw_psi {
                    348:  /arg2 set  %% bft (string)
                    349:  /arg1 set  %% Op  (polynomial)
                    350:  [/bft /bfDt /P /tt /dtt /k /Q /i /m /kk /PPt /PPC /kk /Ss] pushVariables
                    351:  [
                    352:   arg2 /bft set
                    353:   arg1 fw_symbol /P set
                    354:   /bfDt bft xtoDx def
                    355:   /tt bft expand def  /dtt bfDt expand def
                    356:   P bft fw_order /k set
                    357:     << 1 1 k >>
                    358:     {pop tt P mul /P set }
                    359:     for
                    360:     << -1 -1 k >>
                    361:     {pop dtt P mul /P set }
                    362:     for
                    363:   (0) expand /Q set
                    364:   P dtt coefficients 0 get length /m set
                    365:   0 1 << m 1 sub >>
                    366:   {
                    367:     /i set
                    368:     P dtt coefficients 0 get i get /kk set
                    369:     kk (integer) data_conversion /kk set
                    370:     P dtt coefficients 1 get i get /PPt set
                    371:     PPt tt coefficients 1 get 0 get /PPC set
                    372:     bfunction.s expand /Ss set
                    373:     0 1 << kk 1 sub >> {
                    374:       pop
                    375:       PPC Ss mul /PPC set
                    376:       Ss (1) expand sub /Ss set
                    377:     } for
                    378:     Q PPC add /Q set
                    379:   } for
                    380:   Q  /arg1 set
                    381:  ] pop
                    382:  popVariables
                    383:  arg1
                    384: } def
                    385:
                    386: %% get the FW-order
                    387: %% Op (poly) (t) fw_order ---> FW-ord(Op) (integer)
                    388: %% Op should be FW-homogenized.
                    389: /fw_order {
                    390:  /arg2 set  %% bft (string)
                    391:  /arg1 set  %% Op (poly)
                    392:  [/Op /bft /fws /m /fwsDt /k /tt /dtt] pushVariables
                    393:  [
                    394:   arg1 /Op set
                    395:   arg2 /bft set
                    396:   Op fw_symbol /fws set
                    397:   /tt bft expand def
                    398:   /dtt bft xtoDx  expand def
                    399:   fws [[bfunction.s expand  (1).]] replace /fws set
                    400:   fws dtt coefficients 0 get 0 get /m set
                    401:   fws dtt coefficients 1 get 0 get /fwsDt set
                    402:   fwsDt tt coefficients 0 get 0 get /k set
                    403:   m k sub (integer) data_conversion /arg1 set
                    404:  ] pop
                    405:  popVariables
                    406:  arg1
                    407: } def
                    408:
                    409: /remove0 {
                    410:   /arg1 set
                    411:   arg1 (0). eq
                    412:   { } {arg1} ifelse
                    413: } def
                    414:
                    415: %% functions for list operations etc.
                    416:
                    417: /notidentical {
                    418:   /arg2 set
                    419:   /arg1 set
                    420:   arg1 arg2 eq
                    421:   { } {arg1} ifelse
                    422: } def
                    423:
                    424: %% [(x1) (x2) (x3)] ---> (x1,x2,x3)
                    425: /listtostring {
                    426:   /arg1 set
                    427:   [/n /j /ary /str] pushVariables
                    428:   [
                    429:     /ary arg1 def
                    430:     /n ary length def
                    431:     arg1 0 get /str set
                    432:     n 1 gt
                    433:       { str (,) 2 cat_n /str set }{ }
                    434:     ifelse
                    435:     1 1 n 1 sub {
                    436:       /j set
                    437:       j n 1 sub eq
                    438:         {str << ary j get >>  2 cat_n /str set}
                    439:         {str << ary j get >>  (,) 3 cat_n /str set}
                    440:       ifelse
                    441:     } for
                    442:     /arg1 str def
                    443:   ] pop
                    444:   popVariables
                    445:   arg1
                    446: } def
                    447:
                    448: %% (x1) --> (Dx1)
                    449: /xtoDx {
                    450:   /arg1 set
                    451:   @@@.Dsymbol arg1 2 cat_n
                    452: } def
                    453:
                    454: %% concatenate two lists
                    455: /concat {
                    456:   /arg2 set
                    457:   /arg1 set
                    458:   [/n /j /lst1 /lst2 ] pushVariables
                    459:   [
                    460:     /lst1 arg1 def
                    461:     /lst2 arg2 def
                    462:     /n lst2 length def
                    463:     0 1 n 1 sub {
                    464:       /j set
                    465:       lst1 lst2 j get append /lst1 set
                    466:     } for
                    467:     /arg1 lst1 def
                    468:   ] pop
                    469:   popVariables
                    470:   arg1
                    471: } def
                    472:
                    473: %% var (poly) m (integer) ---> var^m (poly)
                    474: /powerZ {
                    475:   /arg2 set %% m
                    476:   /arg1 set %% Var
                    477:   [ /m /var /varstr /pow /nvar] pushVariables
                    478:   [
                    479:     arg1 /var set
                    480:     arg2 /m set
                    481:     var (string) data_conversion /varstr set
                    482:     m -1 gt
                    483:       { var m npower /pow set}
                    484:       { varstr (^(-1)) 2 cat_n expand /nvar set
                    485:         nvar << 0 m sub >> npower /pow set
                    486:        }
                    487:     ifelse
                    488:     pow /arg1 set
                    489:   ] pop
                    490:   popVariables
                    491:   arg1
                    492: } def
                    493:
                    494:
                    495: %% (f) varlist fw_delta ---> [t - f, Dx + f_xDt, ...]
                    496: /fw_delta {
                    497:   /arg2 set %% [(x) (y) ...]
                    498:   /arg1 set %% (f)
                    499:   [ /fstr /f /bft /n /j /varlist /dxvarlist /allvarlist /xi /fxi /dxi /dt
                    500:     /delta /BFdt /BFDtx_weight ] pushVariables
                    501:   [
                    502:     arg1 /fstr set
                    503:     arg2 /varlist set
                    504:     [bfunction.t] varlist join bfunction.a join /allvarlist set
                    505:     bfunction.t xtoDx /BFdt set
                    506:     varlist {xtoDx} map /dxvarlist set
                    507:     varlist length /n set
                    508:     /BFDtx_weight [ [ BFdt 1
                    509:       0 1 n 1 sub {/j set varlist j get 1} for ]
                    510:       [ bfunction.t 1
                    511:       0 1 n 1 sub {/j set dxvarlist j get 1} for ]
                    512:     ] def
                    513:
                    514:     [ allvarlist listtostring ring_of_differential_operators
                    515:      BFDtx_weight weight_vector  0 ] define_ring
                    516:
                    517:     fstr expand /f set
                    518:     bfunction.t expand /bft set
                    519:     BFdt expand /dt set
                    520:     /delta [
                    521:       bft f sub
                    522:       0 1 n 1 sub {
                    523:         /i set
                    524:         varlist i get xtoDx expand /dxi set
                    525:        << dxi f mul >> << f dxi mul >> sub [[(h). (1).]] replace /fxi set
                    526:         dxi << fxi dt mul >> add
                    527:       } for
                    528:     ] def
                    529:     delta {(string) data_conversion} map /arg1 set
                    530:   ] pop
                    531:   popVariables
                    532:   arg1
                    533: } def
                    534:
                    535:
                    536:
                    537:
                    538:

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