Annotation of OpenXM/src/kan96xx/Doc/hol.sm1, Revision 1.2
1.2 ! takayama 1: % $OpenXM$
1.1 maekawa 2: %% hol.sm1, 1998, 11/8, 11/10, 11/14, 11/25, 1999, 5/18, 6/5.
3: %% rank, rrank, characteristic
4: %% This file is error clean.
5: /hol.version (2.990515) def
6: hol.version [(Version)] system_variable gt
7: { [(This package hol.sm1 requires the latest version of kan/sm1) nl
8: (Please get it from http://www.math.kobe-u.ac.jp/KAN)
9: ] cat
10: error
11: } { } ifelse
12:
1.2 ! takayama 13: $hol.sm1, basic package for holonomic systems (C) N.Takayama, 1999, 12/07 $
1.1 maekawa 14: message-quiet
15:
16: /rank.v [(x) (y) (z)] def %% default value of v (variables).
17: /rank.ch [ ] def %% characteristic variety.
18: /rank.verbose 0 def
19: /rank {
20: /arg1 set
21: [/in-rank /aa /typev /setarg /f /v /vsss /vddd
22: /gg /wv /vd /vdweight /chv
23: /one
24: ] pushVariables
25: [(CurrentRingp) (KanGBmessage)] pushEnv
26: [
27:
28: /aa arg1 def
29: aa isArray { } { ( << array >> rank) error } ifelse
30: /setarg 0 def
31: aa { tag } map /typev set
32: typev [ ArrayP ] eq
33: { /f aa 0 get def
34: /v rank.v def
35: /setarg 1 def
36: } { } ifelse
37: typev [ArrayP StringP] eq
38: { /f aa 0 get def
39: /v [ aa 1 get to_records pop ] def
40: /setarg 1 def
41: } { } ifelse
42: typev [ArrayP ArrayP] eq
43: { /f aa 0 get def
44: /v aa 1 get def
45: /setarg 1 def
46: } { } ifelse
47: setarg { } { (rank : Argument mismatch) error } ifelse
48:
49: [(KanGBmessage) rank.verbose ] system_variable
50:
51: f { toString } map /f set
52: v { @@@.Dsymbol 2 1 roll 2 cat_n 1 } map
53: /vddd set %% vddd = [(Dx) 1 (Dy) 1 (Dz) 1]
54: v { @@@.Dsymbol 2 1 roll 2 cat_n } map
55: /vd set %% vd = [(Dx) (Dy) (Dz)]
56: /vdweight
57: vd { [ 2 1 roll -1 ] } map %% vdweight=[[(Dx) -1] [(Dy) -1] [(Dz) -1]]
58: def
59:
60: [v from_records
61: ring_of_differential_operators [vddd] weight_vector 0] define_ring
62: f { . dehomogenize } map /f set
63: [f] groebner_sugar 0 get /gg set
64:
65: /wv vddd weightv def
66: gg { wv init } map /chv set %%obtained the characteristic variety.
67: /rank.ch chv def
68: chv { toString } map /chv set
69:
70: [ v vd join from_records
71: ring_of_polynomials
72: [vddd] vdweight join weight_vector
73: 0
74: ] define_ring
75: [chv {.} map] groebner_sugar 0 get { init } map /chii set
76:
77: /rank.chii chii def
78: rank.verbose { chii message } { } ifelse
79: v {[ 2 1 roll . (1).]} map /one set
80: %% [[(x). (1).] [(y). (1).] [(z). (1).]]
81: %% chii { one replace } map %% buggy code.
82: %% Arg of hilb should be a reduced GB.
83: [chii { one replace } map] groebner 0 get
84: vd hilb /arg1 set
85: ] pop
86: popEnv
87: popVariables
88: arg1
89: } def
90:
91:
92: [(rank)
93: [( a rank b)
94: ( array a; number b)
95: (Example 1 : )
96: $ [ [( (x Dx)^2 + ( y Dy)^2) ( x Dx y Dy -1)] (x,y)] rank :: $
97: (Example 2 : )
98: $[ [( (x^3-y^2) Dx + 3 x^2) ( (x^3-y^2) Dy - 2 y)] (x,y)] rank :: $
99: ]
100: ] putUsages
101: (rank ) messagen-quiet
102:
103: /characteristic.verbose 0 def
104: /characteristic.v [(x) (y) (z)] def
105: /characteristic.ch [ ] def
106: /ch { characteristic } def
107: /characteristic {
108: /arg1 set
109: [/in-rank /aa /typev /setarg /f /v /vsss /vddd
110: /gg /wv /vd /chv
111: /one
112: ] pushVariables
113: [(CurrentRingp) (KanGBmessage)] pushEnv
114: [
115:
116: /aa arg1 def
117: aa isArray { } { ( << array >> characteristic) error } ifelse
118: /setarg 0 def
119: aa { tag } map /typev set
120: typev [ ArrayP ] eq
121: { /f aa 0 get def
122: /v characteristic.v def
123: /setarg 1 def
124: } { } ifelse
125: typev [ArrayP StringP] eq
126: { /f aa 0 get def
127: /v [ aa 1 get to_records pop ] def
128: /setarg 1 def
129: } { } ifelse
130: typev [ArrayP ArrayP] eq
131: { /f aa 0 get def
132: /v aa 1 get def
133: /setarg 1 def
134: } { } ifelse
135: setarg { } { (rank : Argument mismatch) error } ifelse
136:
137: [(KanGBmessage) characteristic.verbose ] system_variable
138:
139: f { toString } map /f set
140: v { @@@.Dsymbol 2 1 roll 2 cat_n 1 } map
141: /vddd set %% vddd = [(Dx) 1 (Dy) 1 (Dz) 1]
142: v { @@@.Dsymbol 2 1 roll 2 cat_n } map
143: /vd set %% vd = [(Dx) (Dy) (Dz)]
144:
145: [v from_records
146: ring_of_differential_operators [vddd] weight_vector 0] define_ring
147: f { . dehomogenize } map /f set
148: [f] groebner_sugar 0 get /gg set
149:
150: /wv vddd weightv def
151: gg { wv init } map /chv set
152: /characteristic.ch [chv] def
153: %% gg { wv init toString} map /chv set %%obtained the characteristic variety.
154: %% /characteristic.ch chv def
155:
156: %% [ v vd join from_records
157: %% ring_of_polynomials
158: %% [vddd] weight_vector
159: %% 0
160: %% ] define_ring
161: %% [chv {.} map] groebner_sugar 0 get /characteristic.ch set
162:
163: characteristic.ch /arg1 set
164: ] pop
165: popEnv
166: popVariables
167: arg1
168: } def
169:
170: [(characteristic)
171: [( a characteristic b)
172: ( array a; number b)
173: (b is the generator of the characteristic variety of a.)
174: (For the algorithm, see Japan J. of Industrial and Applied Math., 1994, 485--497.)
175: (Example 1 : )
176: $ [ [( (x Dx)^2 + ( y Dy)^2) ( x Dx y Dy -1)] (x,y)] characteristic :: $
177: (Example 2 : )
178: $[ [( (x^3-y^2) Dx + 3 x^2) ( (x^3-y^2) Dy - 2 y)] (x,y)] characteristic :: $
179: ]
180: ] putUsages
181: (characteristic ) messagen-quiet
182: [(ch)
183: [(ch is the abbreviation of characteristic.)
184: ( a ch b)
185: ( array a; number b)
186: (b is the generator of the characteristic variety of a.)
187: (For the algorithm, see, Japan J. of Industrial and Applied Math., 1994, 485--497.)
188: (Example 1 : )
189: $ [ [( (x Dx)^2 + ( y Dy)^2) ( x Dx y Dy -1)] (x,y)] ch :: $
190: (Example 2 : )
191: $[ [( (x^3-y^2) Dx + 3 x^2) ( (x^3-y^2) Dy - 2 y)] (x,y)] ch :: $
192: ]
193: ] putUsages
194: (ch ) messagen-quiet
195:
196: %%%% developing rrank.sm1
197: /rrank.v [(x) (y) (z)] def %% default value of v (variables).
198: /rrank.init [ ] def %% initial ideal.
199: /rrank.verbose 0 def
200: /rrank {
201: /arg1 set
202: [/in-rrank /aa /typev /setarg /f /v /vsss /vddd
203: /gg /wv /vd /vdweight
204: /one /i /chv
205: ] pushVariables
206: [(CurrentRingp) (KanGBmessage)] pushEnv
207: [
208:
209: /aa arg1 def
210: aa isArray { } { ( << array >> rrank) error } ifelse
211: /setarg 0 def
212: aa { tag } map /typev set
213: typev [ ArrayP ] eq
214: { /f aa 0 get def
215: /v rrank.v def
216: /setarg 1 def
217: } { } ifelse
218: typev [ArrayP StringP] eq
219: { /f aa 0 get def
220: /v [ aa 1 get to_records pop ] def
221: /setarg 1 def
222: } { } ifelse
223: typev [ArrayP ArrayP] eq
224: { /f aa 0 get def
225: /v aa 1 get def
226: /setarg 1 def
227: } { } ifelse
228: setarg { } { (rrank : Argument mismatch) error } ifelse
229:
230: [(KanGBmessage) rrank.verbose ] system_variable
231:
232: f { toString } map /f set
233: v { @@@.Dsymbol 2 1 roll 2 cat_n 1 } map
234:
235: v { @@@.Dsymbol 2 1 roll 2 cat_n } map
236: /vd set %% vd = [(Dx) (Dy) (Dz)] , v = [(x) (y) (z)]
237: /vdweight
238: [ 0 1 v length 1 sub { /i set v i get << 0 i sub >>
239: vd i get << i >> } for ]
240: def
241: rrank.verbose { vdweight message } { } ifelse
242:
243: [v from_records
244: ring_of_differential_operators [vdweight] weight_vector 0] define_ring
245: f { . dehomogenize homogenize } map /f set
246: [f] groebner 0 get {dehomogenize} map /gg set
247:
248: /wv vdweight weightv def
249: gg { wv init } map /rrank.init set %%obtained the initial ideal
250: rrank.init {toString} map /chv set
251: /arg1 [chv v] rank def
252: ] pop
253: popEnv
254: popVariables
255: arg1
256: } def
257:
258:
259: [(rrank)
260: [( a rrank b)
261: ( array a; number b)
262: (It computes the holonomic rank for regular holonomic system.)
263: (For the algorithm, see Grobner deformations of hypergeometric differential equations, 1999, Springer.)
264: (Chapter 2.)
265: (Example 1 : )
266: $ [ [( (x Dx)^2 + ( y Dy)^2) ( x Dx y Dy -1)] (x,y)] rrank :: $
267: ]
268: ] putUsages
269: (rrank ) messagen-quiet
270:
271: /gb.v 1 def
272: /gb.verbose 0 def
273: /gb {
274: /arg1 set
275: [/in-gb /aa /typev /setarg /f /v
276: /gg /wv /termorder /vec /ans /rr /mm
277: ] pushVariables
278: [(CurrentRingp) (KanGBmessage)] pushEnv
279: [
280:
281: /aa arg1 def
282: aa isArray { } { ( << array >> gb) error } ifelse
283: /setarg 0 def
284: /wv 0 def
285: aa { tag } map /typev set
286: typev [ ArrayP ] eq
287: { /f aa 0 get def
288: /v gb.v def
289: /setarg 1 def
290: } { } ifelse
291: typev [ArrayP StringP] eq
292: { /f aa 0 get def
293: /v aa 1 get def
294: /setarg 1 def
295: } { } ifelse
296: typev [ArrayP ArrayP] eq
297: { /f aa 0 get def
298: /v aa 1 get from_records def
299: /setarg 1 def
300: } { } ifelse
301: typev [ArrayP StringP ArrayP] eq
302: { /f aa 0 get def
303: /v aa 1 get def
304: /wv aa 2 get def
305: /setarg 1 def
306: } { } ifelse
307: typev [ArrayP ArrayP ArrayP] eq
308: { /f aa 0 get def
309: /v aa 1 get from_records def
310: /wv aa 2 get def
311: /setarg 1 def
312: } { } ifelse
313:
314: setarg { } { (gb : Argument mismatch) error } ifelse
315:
316: [(KanGBmessage) gb.verbose ] system_variable
317:
318: %%% Start of the preprocess
319: f getRing /rr set
320: %% To the normal form : matrix expression.
321: f gb.toMatrixOfString /f set
322: /mm gb.itWasMatrix def
323:
324: rr tag 0 eq {
325: %% Define our own ring
326: v isInteger {
327: (Error in gb: Specify variables) error
328: } { } ifelse
329: wv isInteger {
330: [v ring_of_differential_operators
331: 0] define_ring
332: /termorder 1 def
333: }{
334: [v ring_of_differential_operators
335: wv weight_vector
336: 0] define_ring
337: wv gb.isTermOrder /termorder set
338: } ifelse
339: } {
340: %% Use the ring structre given by the input.
341: v isInteger not {
342: (Warning : the given ring definition is not used.) message
343: } { } ifelse
344: rr ring_def
345: /wv rr gb.getWeight def
346: wv gb.isTermOrder /termorder set
347: } ifelse
348: %%% Enf of the preprocess
349:
350:
351: termorder {
352: f { {. dehomogenize} map } map /f set
353: [f] groebner_sugar 0 get /gg set
354: }{
355: f { {. dehomogenize} map} map /f set
356: f fromVectors { homogenize } map /f set
357: [f] groebner 0 get /gg set
358: }ifelse
359: wv isInteger {
360: /ans [gg gg {init} map] def
361: }{
362: /ans [gg gg {wv 0 get weightv init} map] def
363: }ifelse
364:
365: %% Postprocess : recover the matrix expression.
366: mm {
367: ans { /tmp set [mm tmp] toVectors } map
368: /ans set
369: }{ }
370: ifelse
371: %%
372:
373: /arg1 ans def
374: ] pop
375: popEnv
376: popVariables
377: arg1
378: } def
379: (gb ) messagen-quiet
380:
381: /pgb {
382: /arg1 set
383: [/in-pgb /aa /typev /setarg /f /v
384: /gg /wv /termorder /vec /ans /rr /mm
385: ] pushVariables
386: [(CurrentRingp) (KanGBmessage) (UseCriterion1)] pushEnv
387: [
388:
389: /aa arg1 def
390: aa isArray { } { (<< array >> pgb) error } ifelse
391: /setarg 0 def
392: /wv 0 def
393: aa { tag } map /typev set
394: typev [ ArrayP ] eq
395: { /f aa 0 get def
396: /v gb.v def
397: /setarg 1 def
398: } { } ifelse
399: typev [ArrayP StringP] eq
400: { /f aa 0 get def
401: /v aa 1 get def
402: /setarg 1 def
403: } { } ifelse
404: typev [ArrayP ArrayP] eq
405: { /f aa 0 get def
406: /v aa 1 get from_records def
407: /setarg 1 def
408: } { } ifelse
409: typev [ArrayP StringP ArrayP] eq
410: { /f aa 0 get def
411: /v aa 1 get def
412: /wv aa 2 get def
413: /setarg 1 def
414: } { } ifelse
415: typev [ArrayP ArrayP ArrayP] eq
416: { /f aa 0 get def
417: /v aa 1 get from_records def
418: /wv aa 2 get def
419: /setarg 1 def
420: } { } ifelse
421:
422: setarg { } { (pgb : Argument mismatch) error } ifelse
423:
424: [(KanGBmessage) gb.verbose ] system_variable
425:
426: %%% Start of the preprocess
427: f getRing /rr set
428: %% To the normal form : matrix expression.
429: f gb.toMatrixOfString /f set
430: /mm gb.itWasMatrix def
431:
432: rr tag 0 eq {
433: %% Define our own ring
434: v isInteger {
435: (Error in pgb: Specify variables) error
436: } { } ifelse
437: wv isInteger {
438: [v ring_of_polynomials
439: 0] define_ring
440: /termorder 1 def
441: }{
442: [v ring_of_polynomials
443: wv weight_vector
444: 0] define_ring
445: wv gb.isTermOrder /termorder set
446: } ifelse
447: } {
448: %% Use the ring structre given by the input.
449: v isInteger not {
450: (Warning : the given ring definition is not used.) message
451: } { } ifelse
452: rr ring_def
453: /wv rr gb.getWeight def
454: wv gb.isTermOrder /termorder set
455: } ifelse
456: %%% Enf of the preprocess
457:
458:
459: termorder {
460: f { {. dehomogenize} map } map /f set
461: [(UseCriterion1) 1] system_variable
462: [f] groebner_sugar 0 get /gg set
463: [(UseCriterion1) 0] system_variable
464: }{
465: f { {. dehomogenize} map} map /f set
466: f fromVectors { homogenize } map /f set
467: [(UseCriterion1) 1] system_variable
468: [f] groebner 0 get /gg set
469: [(UseCriterion1) 0] system_variable
470: }ifelse
471: wv isInteger {
472: /ans [gg gg {init} map] def
473: }{
474: /ans [gg gg {wv 0 get weightv init} map] def
475: }ifelse
476:
477: %% Postprocess : recover the matrix expression.
478: mm {
479: ans { /tmp set [mm tmp] toVectors } map
480: /ans set
481: }{ }
482: ifelse
483: %%
484:
485: /arg1 ans def
486: ] pop
487: popEnv
488: popVariables
489: arg1
490: } def
491:
492: /pgb.old {
493: /arg1 set
494: [/in-pgb /aa /typev /setarg /f /v
495: /gg /wv /termorder /vec /ans
496: ] pushVariables
497: [(CurrentRingp) (KanGBmessage) (UseCriterion1)] pushEnv
498: [
499:
500: /aa arg1 def
501: aa isArray { } { (array pgb) message (pgb) usage error } ifelse
502: /setarg 0 def
503: /wv 0 def
504: aa { tag } map /typev set
505: typev [ ArrayP ] eq
506: { /f aa 0 get def
507: /v gb.v def
508: /setarg 1 def
509: } { } ifelse
510: typev [ArrayP StringP] eq
511: { /f aa 0 get def
512: /v aa 1 get def
513: /setarg 1 def
514: } { } ifelse
515: typev [ArrayP ArrayP] eq
516: { /f aa 0 get def
517: /v aa 1 get from_records def
518: /setarg 1 def
519: } { } ifelse
520: typev [ArrayP StringP ArrayP] eq
521: { /f aa 0 get def
522: /v aa 1 get def
523: /wv aa 2 get def
524: /setarg 1 def
525: } { } ifelse
526: typev [ArrayP ArrayP ArrayP] eq
527: { /f aa 0 get def
528: /v aa 1 get from_records def
529: /wv aa 2 get def
530: /setarg 1 def
531: } { } ifelse
532:
533: setarg { } { (pgb : Argument mismatch) message error } ifelse
534:
535: [(KanGBmessage) gb.verbose ] system_variable
536:
537: %% Input must not be vectors.
538: f { toString } map /f set
539:
540: wv isInteger {
541: [v ring_of_polynomials
542: 0] define_ring
543: /termorder 1 def
544: }{
545: [v ring_of_polynomials
546: wv weight_vector
547: 0] define_ring
548: wv gb.isTermOrder /termorder set
549: } ifelse
550: termorder {
551: f { . dehomogenize } map /f set
552: [(UseCriterion1) 1] system_variable
553: [f] groebner_sugar 0 get /gg set
554: [(UseCriterion1) 0] system_variable
555: }{
556: f { . dehomogenize homogenize} map /f set
557: [(UseCriterion1) 1] system_variable
558: [f] groebner 0 get /gg set
559: [(UseCriterion1) 0] system_variable
560: }ifelse
561: wv isInteger {
562: /ans [gg gg {init} map] def
563: }{
564: /ans [gg gg {wv 0 get weightv init} map] def
565: }ifelse
566: /arg1 ans def
567: ] pop
568: popEnv
569: popVariables
570: arg1
571: } def
572: (pgb ) messagen-quiet
573:
574: /gb.toMatrixOfString {
575: /arg1 set
576: [/in-gb.toMatrixOfString /ff /aa /ans] pushVariables
577: [
578: /aa arg1 def
579: aa length 0 eq { /ans [ ] def /gb.toMatrixOfString.LLL goto }{ } ifelse
580: aa 0 get isArray {
581: /gb.itWasMatrix aa 0 get length def
582: }{
583: /gb.itWasMatrix 0 def
584: } ifelse
585: aa {
586: /ff set
587: ff isArray {
588: ff {toString} map /ff set
589: }{
590: [ff toString] /ff set
591: } ifelse
592: ff
593: } map /ans set
594: /gb.toMatrixOfString.LLL
595: /arg1 ans def
596: ] pop
597: popVariables
598: arg1
599: } def
600: [(gb.toMatrixOfString)
601: [(It translates given input into a matrix form which is a data structure)
602: (for computations of kernel, image, cokernel, etc.)
603: (gb.itWasMatrix is set to the length of the input vector.)
604: $Example 1: $
605: $ [ (x). (y).] gb.toMatrixOfString ==> [[(x)] [(y)]] $
606: $ gb.itWasMatrix is 0.$
607: $Example 2: $
608: $ [ [(x). (1).] [(y). (0).]] gb.toMatrixOfString ==> [ [(x) (1)] [(y) (0)]] $
609: $ gb.itWasMatrix is 2.$
610: ]] putUsages
611:
612: /gb.toMatrixOfPoly {
613: /arg1 set
614: [/in-gb.toMatrixOfPoly /ff /aa /ans] pushVariables
615: [
616: /aa arg1 def
617: aa length 0 eq { /ans [ ] def /gb.toMatrixOfPoly.LLL goto }{ } ifelse
618: aa 0 get isArray {
619: /gb.itWasMatrix aa 0 get length def
620: }{
621: /gb.itWasMatrix 0 def
622: } ifelse
623: aa {
624: /ff set
625: ff isArray {
626: }{
627: [ff] /ff set
628: } ifelse
629: ff
630: } map /ans set
631: /gb.toMatrixOfPoly.LLL
632: /arg1 ans def
633: ] pop
634: popVariables
635: arg1
636: } def
637: [(gb.toMatrixOfPoly)
638: [(It translates given input into a matrix form which is a data structure)
639: (for computations of kernel, image, cokernel, etc.)
640: (gb.itWasMatrix is set to the length of the input vector.)
641: $Example 1: $
642: $ [ (x). (y).] gb.toMatrixOfPoly ==> [[(x)] [(y)]] $
643: $ gb.itWasMatrix is 0.$
644: $Example 2: $
645: $ [ [(x). (1).] [(y). (0).]] gb.toMatrixOfPoly ==> [ [(x) (1)] [(y) (0)]] $
646: $ gb.itWasMatrix is 2.$
647: ]] putUsages
648:
649: /gb.getWeight {
650: /arg1 set
651: [/in-gb.getWeight /rr /ww /vv /ans /nn /ii] pushVariables
652: [(CurrentRingp)] pushEnv
653: [
654: /rr arg1 def
655: rr ring_def
656: getVariableNames /vv set
657: [(orderMatrix)] system_variable 0 get /ww set
658: /nn vv length 1 sub def
659: [0 1 nn {
660: /ii set
661: ww ii get 0 eq {
662: } {
663: vv ii get
664: ww ii get
665: } ifelse
666: } for
667: ] /ans set
668: /arg1 [ans] def
669: ] pop
670: popEnv
671: popVariables
672: arg1
673: } def
674: [(gb.getWeight)
675: [(ring gb.getWeight wv)
676: (It gets the weight vector field of the ring ring.)
677: ]] putUsages
678:
679:
680: /gb.isTermOrder {
681: /arg1 set
682: [/in-gb.isTermOrder /vv /ww /yes /i /j] pushVariables
683: [
684: /vv arg1 def
685: /yes 1 def
686: 0 1 vv length 1 sub {
687: /i set
688: /ww vv i get def
689: 0 1 ww length 1 sub {
690: /j set
691: ww j get isInteger {
692: ww j get 0 lt { /yes 0 def } { } ifelse
693: }{ } ifelse
694: }for
695: }for
696: /arg1 yes def
697: ] pop
698: popVariables
699: arg1
700: } def
701: [(gb)
702: [(a gb b)
703: (array a; array b;)
704: (b : [g ii]; array g; array in; g is a Grobner basis of f)
705: ( in the ring of differential operators.)
706: $ ii is the initial ideal in case of w is given or <<a>> belongs$
707: $ to a ring. In the other cases, it returns the initial monominal.$
708: (a : [f ]; array f; f is a set of generators of an ideal in a ring.)
709: (a : [f v]; array f; string v; v is the variables. )
710: (a : [f v w]; array f; string v; array of array w; w is the weight matirx.)
711: ( )
712: $Example 1: [ [( (x Dx)^2 + (y Dy)^2 -1) ( x y Dx Dy -1)] (x,y) $
713: $ [ [ (Dx) 1 ] ] ] gb pmat ; $
714: (Example 2: )
715: (To put h=1, type in, e.g., )
716: $ [ [(2 x Dx + 3 y Dy+6) (2 y Dx + 3 x^2 Dy)] (x,y) $
717: $ [[(x) -1 (Dx) 1 (y) -1 (Dy) 1]]] gb /gg set gg dehomogenize pmat ;$
718: ( )
719: $Example 3: [ [( (x Dx)^2 + (y Dy)^2 -1) ( x y Dx Dy -1)] (x,y) $
720: $ [ [ (Dx) 1 (Dy) 1] ] ] gb pmat ; $
721: ( )
722: $Example 4: [[ [(x^2) (y+x)] [(x+y) (y^3)] [(2 x^2+x y) (y+x+x y^3)]] (x,y) $
723: $ [ [ (x) -1 (y) -1] ] ] gb pmat ; $
724: ( )
725: (cf. gb, groebner, groebner_sugar, syz. )
726: ]] putUsages
727:
728: [(pgb)
729: [(a pgb b)
730: (array a; array b;)
731: (b : [g ii]; array g; array in; g is a Grobner basis of f)
732: ( in the ring of polynomials.)
733: $ ii is the initial ideal in case of w is given or <<a>>belongs$
734: $ to a ring. In the other cases, it returns the initial monominal.$
735: (a : [f ]; array f; f is a set of generators of an ideal in a ring.)
736: (a : [f v]; array f; string v; v is the variables.)
737: (a : [f v w]; array f; string v; array of array w; w is the weight matirx.)
738: $Example 1: [(x,y) ring_of_polynomials 0] define_ring $
739: $ [ [(x^2+y^2-4). (x y -1).] ] pgb :: $
740: $Example 2: [ [(x^2+y^2) (x y)] (x,y) [ [(x) -1 (y) -1] ] ] pgb :: $
741: (cf. gb, groebner, groebner_sugar, syz. )
742: ]] putUsages
743:
744:
745: %/syz.v 1 def
746: /syz.v 1 def
747: /syz.verbose 0 def
748: /syz {
749: /arg1 set
750: [/in-syz /aa /typev /setarg /f /v
751: /gg /wv /termorder /vec /ans /ggall /vectorInput /vsize /gtmp /gtmp2
752: /rr /mm
753: ] pushVariables
754: [(CurrentRingp) (KanGBmessage)] pushEnv
755: [
756:
757: /aa arg1 def
758: aa isArray { } { (<< array >> syz) error } ifelse
759: /setarg 0 def
760: /wv 0 def
761: aa { tag } map /typev set
762: typev [ ArrayP ] eq
763: { /f aa 0 get def
764: /v syz.v def
765: /setarg 1 def
766: } { } ifelse
767: typev [ArrayP StringP] eq
768: { /f aa 0 get def
769: /v aa 1 get def
770: /setarg 1 def
771: } { } ifelse
772: typev [ArrayP ArrayP] eq
773: { /f aa 0 get def
774: /v aa 1 get from_records def
775: /setarg 1 def
776: } { } ifelse
777: typev [ArrayP StringP ArrayP] eq
778: { /f aa 0 get def
779: /v aa 1 get def
780: /wv aa 2 get def
781: /setarg 1 def
782: } { } ifelse
783: typev [ArrayP ArrayP ArrayP] eq
784: { /f aa 0 get def
785: /v aa 1 get from_records def
786: /wv aa 2 get def
787: /setarg 1 def
788: } { } ifelse
789:
790: setarg { } { (syz : Argument mismatch) error } ifelse
791:
792: [(KanGBmessage) syz.verbose ] system_variable
793:
794:
795:
796: %%% Start of the preprocess
797: f getRing /rr set
798: %% To the normal form : matrix expression.
799: f gb.toMatrixOfString /f set
800: /mm gb.itWasMatrix def
801: mm 0 gt {
802: /vectorInput 1 def
803: }{
804: /vectorInput 1 def
805: } ifelse
806:
807: rr tag 0 eq {
808: %% Define our own ring
809: v isInteger {
810: (Error in syz: Specify variables) error
811: } { } ifelse
812: wv isInteger {
813: [v ring_of_differential_operators
814: 0] define_ring
815: /termorder 1 def
816: }{
817: [v ring_of_differential_operators
818: wv weight_vector
819: 0] define_ring
820: wv gb.isTermOrder /termorder set
821: } ifelse
822: }{
823: %% Use the ring structre given by the input.
824: v isInteger not {
825: (Warning : the given ring definition is not used.) message
826: } { } ifelse
827: rr ring_def
828: /wv rr gb.getWeight def
829: wv gb.isTermOrder /termorder set
830: } ifelse
831: %%% Enf of the preprocess
832:
833: termorder {
834: f { {. dehomogenize} map } map /f set
835: [f [(needBack) (needSyz)]] groebner_sugar /ggall set
836: ggall 2 get /gg set
837: }{
838: f { {. dehomogenize } map homogenize } map /f set
839: [f [(needBack) (needSyz)]] groebner /ggall set
840: ggall 2 get /gg set
841: }ifelse
842: vectorInput {
843: /vsize f 0 get length def %% input vector size.
844: /gtmp ggall 0 get def
845: [vsize gtmp] toVectors /gtmp set
846: ggall 0 gtmp put
847: }{ } ifelse
848: /arg1 [gg dehomogenize ggall] def
849: ] pop
850: popEnv
851: popVariables
852: arg1
853: } def
854: (syz ) messagen-quiet
855:
856: [(syz)
857: [(a syz [b c])
858: (array a; array b; array c)
859: (b is a set of generators of the syzygies of f.)
860: (c = [gb, backward transformation, syzygy without dehomogenization].)
861: (See groebner.)
862: (a : [f ]; array f; f is a set of generators of an ideal in a ring.)
863: (a : [f v]; array f; string v; v is the variables.)
864: (a : [f v w]; array f; string v; array of array w; w is the weight matirx.)
865: $Example 1: [(x,y) ring_of_polynomials 0] define_ring $
866: $ [ [(x^2+y^2-4). (x y -1).] ] syz :: $
867: $Example 2: [ [(x^2+y^2) (x y)] (x,y) [ [(x) -1 (y) -1] ] ] syz :: $
868: $Example 3: [ [( (x Dx)^2 + (y Dy)^2 -1) ( x y Dx Dy -1)] (x,y) $
869: $ [ [ (Dx) 1 ] ] ] syz pmat ; $
870: $Example 4: [ [(2 x Dx + 3 y Dy+6) (2 y Dx + 3 x^2 Dy)] (x,y) $
871: $ [[(x) -1 (Dx) 1 (y) -1 (Dy) 1]]] syz pmat ;$
872: $Example 5: [ [ [(x^2) (y+x)] [(x+y) (y^3)] [(2 x^2+x y) (y+x+x y^3)]] $
873: $ (x,y) ] syz pmat ;$
874: $Example 6: [ [ [(x^2) (y+x)] [(x+y) (y^3)] [(2 x^2+x y) (y+x+x y^3)]] $
875: $ (x,y) [[(x) -1 (y) -2]] ] syz pmat ;$
876: $Example 7: [ [ [(0) (0)] [(0) (0)] [(x) (y)]] $
877: $ [(x) (y)]] syz pmat ;$
878: ]] putUsages
879:
880:
881: %%%%%%%%%%%%%%%%%% package fs %%%%%%%%%%%%%%%%%%%%%%%
882: [(genericAnn)
883: [ (f [s v1 v2 ... vn] genericAnn [L1 ... Lm])
884: (L1, ..., Lm are annihilating ideal for f^s.)
885: (f is a polynomial of v1, ..., vn)
886: (<string> | <poly> f, s, v1, ..., vn ; <poly> L1, ..., Lm )
887: $Example: (x^3+y^3+z^3) [(s) (x) (y) (z)] genericAnn$
888: ]
889: ] putUsages ( genericAnn ) messagen-quiet
890: /fs.verbose 0 def
891: /genericAnn {
892: /arg2 set /arg1 set
893: [/in-genericAnn /f /vlist /s /vvv /nnn /rrr
894: /v1 /ops /ggg /ggg0
895: ] pushVariables
896: [(CurrentRingp) (KanGBmessage)] pushEnv
897: [
898: /f arg1 def /vlist arg2 def
899: f toString /f set
900: vlist { toString } map /vlist set
901: [(KanGBmessage) fs.verbose] system_variable
902: /s vlist 0 get def
903: /vvv (_u,_v,_t,) vlist rest { (,) 2 cat_n } map aload length /nnn set
904: s nnn 2 add cat_n def
905: fs.verbose { vvv message } { }ifelse
906: [vvv ring_of_differential_operators
907: [[(_u) 1 (_v) 1]] weight_vector 0] define_ring /rrr set
908:
909: [ (_u*_t). f . sub (_u*_v-1). ]
910: vlist rest { /v1 set
911: %%D-clean f . (D) v1 2 cat_n . 1 diff0 (_v*D_t). mul
912: f . @@@.Dsymbol v1 2 cat_n . 1 diff0 [(_v*) @@@.Dsymbol (_t)] cat . mul
913: @@@.Dsymbol v1 2 cat_n . add } map
914: join
915: /ops set
916: ops {[[(h). (1).]] replace } map /ops set
917: fs.verbose { ops message } { }ifelse
918: [ops] groebner_sugar 0 get /ggg0 set
919: fs.verbose { ggg0 message } { } ifelse
920: ggg0 [(_u) (_v)] eliminatev
921: %%D-clean { [(_t).] [ (D_t).] [s .] distraction
922: { [(_t).] [ [@@@.Dsymbol (_t)] cat .] [s .] distraction
923: [[s . << (0). s . sub (1). sub >>]] replace
924: } map /arg1 set
925: ] pop
926: popEnv
927: popVariables
928: arg1
929: } def
930:
931: %% Find differential equations for f^(m), r0 the minimal integral root.
932: [(annfs)
933: [( [ f v m r0] annfs g )
934: (It returns the annihilating ideal of f^m where r0 must be smaller)
935: (or equal to the minimal integral root of the b-function.)
936: (Or, it returns the annihilating ideal of f^r0, r0 and the b-function)
937: (where r0 is the minial integral root of b.)
938: (For the algorithm, see J. Pure and Applied Algebra 117&118(1997), 495--518.)
939: (Example 1: [(x^2+y^2+z^2+t^2) (x,y,z,t) -1 -2] annfs :: )
940: $ It returns the annihilating ideal of (x^2+y^2+z^2+t^2)^(-1).$
941: (Example 2: [(x^2+y^2+z^2+t^2) (x,y,z,t)] annfs :: )
942: $ It returns the annihilating ideal of f^r0 and [r0, b-function]$
943: $ where r0 is the minimal integral root of the b-function.$
944: (Example 3: [(x^2+y^2+z^2) (x,y,z) -1 -1] annfs :: )
945: (Example 4: [(x^3+y^3+z^3) (x,y,z)] annfs :: )
946: (Example 5: [((x1+x2+x3)(x1 x2 + x2 x3 + x1 x3) - t x1 x2 x3 ) )
947: ( (t,x1,x2,x3) -1 -2] annfs :: )
948: ( Note that the example 4 uses huge memory space.)
949: ]] putUsages
950: ( annfs ) messagen-quiet
951: /annfs.verbose fs.verbose def
952: /annfs.v [(x) (y) (z)] def
953: /annfs.s (_s) def
954: %% The first variable must be s.
955: /annfs {
956: /arg1 set
957: [/in-annfs /aa /typev /setarg /v /m /r0 /gg /ss /fs /gg2
958: /ans /vtmp /w2 /velim /bbb /rrr /r0
959: ] pushVariables
960: [(CurrentRingp) (KanGBmessage)] pushEnv
961: [
962:
963: /aa arg1 def
964: aa isArray { } { ( << array >> annfs) error } ifelse
965: /setarg 0 def
966: aa { tag } map /typev set
967: /r0 [ ] def
968: /m [ ] def
969: /v annfs.v def
970: aa 0 << aa 0 get toString >> put
971: typev [ StringP ] eq
972: { /f aa 0 get def
973: /setarg 1 def
974: } { } ifelse
975: typev [StringP StringP] eq
976: { /f aa 0 get def
977: /v [ aa 1 get to_records pop ] def
978: /setarg 1 def
979: } { } ifelse
980: typev [StringP ArrayP] eq
981: { /f aa 0 get def
982: /v aa 1 get def
983: /setarg 1 def
984: } { } ifelse
985: typev [StringP ArrayP IntegerP IntegerP] eq
986: { /f aa 0 get def
987: /v aa 1 get def
988: /m aa 2 get def
989: /r0 aa 3 get def
990: /setarg 1 def
991: } { } ifelse
992: typev [StringP StringP IntegerP IntegerP] eq
993: { /f aa 0 get def
994: /v [ aa 1 get to_records pop ] def
995: /m aa 2 get def
996: /r0 aa 3 get def
997: /setarg 1 def
998: } { } ifelse
999: setarg 1 eq { } { (annfs : wrong argument) error } ifelse
1000:
1001: [annfs.s] v join /v set
1002:
1003: /ss v 0 get def
1004: annfs.verbose {
1005: (f, v, s, f^{m}, m+r0 = ) messagen
1006: [ f (, ) v (, ) ss (, )
1007: (f^) m (,) m (+) r0 ] {messagen} map ( ) message
1008: } { } ifelse
1009:
1010: f v genericAnn /fs set
1011:
1012: annfs.verbose {
1013: (genericAnn is ) messagen fs message
1014: } { } ifelse
1015: [(KanGBmessage) annfs.verbose] system_variable
1016:
1017: m isArray {
1018: %% Now, let us find the b-function. /vtmp /w2 /velim /bbb /rrr /r0
1019: v rest { /vtmp set vtmp @@@.Dsymbol vtmp 2 cat_n } map /velim set
1020: velim { 1 } map /w2 set
1021: annfs.verbose { w2 message } { } ifelse
1022: [v from_records ring_of_differential_operators
1023: [w2] weight_vector 0] define_ring
1024: [ fs { toString . } map [ f toString . ] join ]
1025: groebner_sugar 0 get velim eliminatev 0 get /bbb set
1026: [[(s) annfs.s] from_records ring_of_polynomials 0] define_ring
1027: bbb toString . [[annfs.s . (s).]] replace /bbb set
1028: annfs.verbose { bbb message } { } ifelse
1029: bbb findIntegralRoots /rrr set
1030: rrr 0 get /r0 set %% minimal integral root.
1031: annfs.verbose { rrr message } { } ifelse
1032: fs 0 get (ring) dc ring_def
1033: fs { [[annfs.s . r0 toString .]] replace } map /ans set
1034: /ans [ans [r0 bbb]] def
1035: /annfs.label1 goto
1036: } { } ifelse
1037: m 0 ge {
1038: (annfs works only for getting annihilating ideal for f^(negative))
1039: error
1040: } { } ifelse
1041: r0 isArray {
1042: [(Need to compute the minimal root of b-function) nl
1043: (It has not been implemented.) ] cat
1044: error
1045: } { } ifelse
1046:
1047: [v from_records ring_of_differential_operators 0] define_ring
1048: fs {toString . dehomogenize [[ss . r0 (poly) dc]] replace}
1049: map /gg set
1050: annfs.verbose { gg message } { } ifelse
1051:
1052: [ [f . << m r0 sub >> npower ] gg join
1053: [(needBack) (needSyz)]] groebner_sugar 2 get /gg2 set
1054:
1055: gg2 { 0 get } map /ans set
1056: /ans ans { dup (0). eq {pop} { } ifelse } map def
1057:
1058: /annfs.label1
1059: /arg1 ans def
1060: ] pop
1061: popEnv
1062: popVariables
1063: arg1
1064: } def
1065:
1066: /genericAnnWithL.s (s) def
1067: /annfs.verify 0 def
1068: /genericAnnWithL {
1069: /arg1 set
1070: [/in-genericAnnWithL /aa /typev /setarg /v /m /r0 /gg /ss /fs /gg2
1071: /ans /vtmp /w2 /velim /bbb /rrr /r0 /myL /mygb /jj
1072: ] pushVariables
1073: [(CurrentRingp) (KanGBmessage) (Homogenize)] pushEnv
1074: [
1075:
1076: /aa arg1 def
1077: aa isArray { } { ( << array >> annfs) error } ifelse
1078: /setarg 0 def
1079: aa { tag } map /typev set
1080: /r0 [ ] def
1081: /m [ ] def
1082: /v annfs.v def
1083: aa 0 << aa 0 get toString >> put
1084: typev [ StringP ] eq
1085: { /f aa 0 get def
1086: /setarg 1 def
1087: } { } ifelse
1088: typev [StringP StringP] eq
1089: { /f aa 0 get def
1090: /v [ aa 1 get to_records pop ] def
1091: /setarg 1 def
1092: } { } ifelse
1093: typev [StringP ArrayP] eq
1094: { /f aa 0 get def
1095: /v aa 1 get def
1096: /setarg 1 def
1097: } { } ifelse
1098: setarg 1 eq { } { (genericAnnWithL : wrong argument) error } ifelse
1099:
1100: [genericAnnWithL.s] v join /v set
1101:
1102: /ss v 0 get def
1103: annfs.verbose {
1104: (f, v, s, f^{m}, m+r0 = ) messagen
1105: [ f (, ) v (, ) ss (, )
1106: (f^) m (,) m (+) r0 ] {messagen} map ( ) message
1107: } { } ifelse
1108:
1109: f v genericAnn /fs set
1110:
1111: annfs.verbose {
1112: (genericAnn is ) messagen fs message
1113: } { } ifelse
1114: [(KanGBmessage) annfs.verbose] system_variable
1115:
1116: m isArray {
1117: %% Now, let us find the b-function. /vtmp /w2 /velim /bbb /rrr /r0
1118: v rest { /vtmp set vtmp @@@.Dsymbol vtmp 2 cat_n } map /velim set
1119: velim { 1 } map /w2 set
1120: annfs.verbose { w2 message } { } ifelse
1121: [v from_records ring_of_differential_operators
1122: [w2] weight_vector 0] define_ring
1123:
1124: [ [ f toString . ] fs { toString . } map join [(needBack)]]
1125: groebner_sugar /mygb set
1126: mygb 0 get velim eliminatev 0 get /bbb set
1127: mygb 0 get bbb position /jj set
1128: mygb 1 get jj get 0 get /myL set
1129:
1130: annfs.verbose { bbb message } { } ifelse
1131:
1132: annfs.verify {
1133: (Verifying L f - b belongs to genericAnn(f)) message
1134: [(Homogenize) 0] system_variable
1135: << myL f . mul bbb sub >>
1136: [fs { toString . } map] groebner_sugar 0 get
1137: reduction 0 get message
1138: (Is it zero? Then it's fine.) message
1139: } { } ifelse
1140:
1141: /ans [bbb [myL fs] ] def
1142: /annfs.label1 goto
1143: } { } ifelse
1144:
1145: /annfs.label1
1146: /arg1 ans def
1147: ] pop
1148: popEnv
1149: popVariables
1150: arg1
1151: } def
1152:
1153:
1154: [(genericAnnWithL)
1155: [$[f v] genericAnnWithL [b [L I]]$
1156: $String f,v; poly b,L; array of poly I;$
1157: $f is a polynomial given by a string. v is the variables.$
1158: $ v must not contain names s, e.$
1159: $b is the b-function (Bernstein-Sato polynomial) for f and$
1160: $ L is the operator satisfying L f^{s+1} = b(s) f^s $
1161: $ I is the annihilating ideal of f^s.$
1162: $cf. bfunction, annfs, genericAnn.$
1163: $Example 1: [(x^2+y^2) (x,y)] genericAnnWithL ::$
1164: $Example 2: [(x^2+y^2+z^2) (x,y,z)] genericAnnWithL ::$
1165: $Example 3: [(x^3-y^2 z^2) (x,y,z)] genericAnnWithL ::$
1166: ]] putUsages
1.2 ! takayama 1167:
! 1168: /reduction*.noH 0 def
! 1169: /reduction* {
! 1170: /arg1 set
! 1171: [/in-reduction* /aa /typev /setarg /f /v
! 1172: /gg /wv /termorder /vec /ans /rr /mm /h /size /a0 /a3
! 1173: ] pushVariables
! 1174: [(CurrentRingp) (KanGBmessage)] pushEnv
! 1175: [
! 1176:
! 1177: /aa arg1 def
! 1178: aa isArray { } { ( << array >> reduction*) error } ifelse
! 1179: /setarg 0 def
! 1180: /wv 0 def
! 1181: aa { tag } map /typev set
! 1182: typev [StringP ArrayP ArrayP] eq
! 1183: typev [ArrayP ArrayP ArrayP] eq or
! 1184: typev [PolyP ArrayP ArrayP] eq or
! 1185: { /h aa 0 get def
! 1186: /f aa 1 get def
! 1187: /v aa 2 get from_records def
! 1188: /setarg 1 def
! 1189: } { } ifelse
! 1190: typev [StringP ArrayP ArrayP ArrayP] eq
! 1191: typev [ArrayP ArrayP ArrayP ArrayP] eq or
! 1192: typev [PolyP ArrayP ArrayP ArrayP] eq or
! 1193: { /h aa 0 get def
! 1194: /f aa 1 get def
! 1195: /v aa 2 get from_records def
! 1196: /wv aa 3 get def
! 1197: /setarg 1 def
! 1198: } { } ifelse
! 1199:
! 1200: setarg { } { (reduction* : Argument mismatch) error } ifelse
! 1201:
! 1202: [(KanGBmessage) gb.verbose ] system_variable
! 1203:
! 1204: %%% Start of the preprocess
! 1205: f getRing /rr set
! 1206:
! 1207:
! 1208: rr tag 0 eq {
! 1209: %% Define our own ring
! 1210: v isInteger {
! 1211: (Error in reduction*: Specify variables) error
! 1212: } { } ifelse
! 1213: wv isInteger {
! 1214: [v ring_of_differential_operators
! 1215: 0] define_ring
! 1216: /termorder 1 def
! 1217: }{
! 1218: [v ring_of_differential_operators
! 1219: wv weight_vector
! 1220: 0] define_ring
! 1221: wv gb.isTermOrder /termorder set
! 1222: } ifelse
! 1223: } {
! 1224: %% Use the ring structre given by the input.
! 1225: v isInteger not {
! 1226: (Warning : the given ring definition is not used.) message
! 1227: } { } ifelse
! 1228: rr ring_def
! 1229: /wv rr gb.getWeight def
! 1230: wv gb.isTermOrder /termorder set
! 1231: } ifelse
! 1232: %%% Enf of the preprocess
! 1233:
! 1234: f 0 get isArray {
! 1235: /size f 0 get length def
! 1236: f { { toString . } map } map /f set
! 1237: f fromVectors /f set
! 1238: }{
! 1239: /size -1 def
! 1240: f { toString . } map /f set
! 1241: } ifelse
! 1242:
! 1243: h isArray {
! 1244: h { toString . } map /h set
! 1245: [h] fromVectors 0 get /h set
! 1246: }{
! 1247: h toString . /h set
! 1248: } ifelse
! 1249: f { toString . } map /f set
! 1250: reduction*.noH {
! 1251: h f reduction-noH /ans set
! 1252: } {
! 1253: h f reduction /ans set
! 1254: } ifelse
! 1255: size -1 eq not {
! 1256: [size ans 0 get] toVectors /a0 set
! 1257: [size ans 3 get] toVectors /a3 set
! 1258: /ans [a0 ans 1 get ans 2 get a3] def
! 1259: } { } ifelse
! 1260: /arg1 ans def
! 1261: ] pop
! 1262: popEnv
! 1263: popVariables
! 1264: arg1
! 1265: } def
! 1266:
! 1267:
! 1268: [(reduction*)
! 1269: [([f base v] reduction* [h c0 syz input])
! 1270: ([f base v weight] reduction* [h c0 syz input])
! 1271: (reduction* is an user interface for reduction and reduction-noH.)
! 1272: (If reduction*.noH is one, then reduction-noH will be called.)
! 1273: (Example 1: [(x^2) [(x^2+y^2-4) (x y-1)] [(x) (y)]] reduction* )
! 1274: (Example 2: [[(1) (y^2-1)] [ [(0) (y-1)] [(1) (y+1)]] [(x) (y)]] reduction*)
! 1275: (Example 3: [(x^2) [(x^2+y^2-4) (x y-1)] [(x) (y)] [[(x) 10]] ] reduction* )
! 1276: ]] putUsages
1.1 maekawa 1277:
1278: ( ) message-quiet ;
1279:
1280:
1281:
1282:
1283:
1284:
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to hol.sm1 CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM / src / kan96xx / Doc