![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to generic_polynomials.adb CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / PHC / Ada / Math_Lib / Polynomials|
Return to generic_polynomials.adb CVS log | Up to [local] / OpenXM_contrib / PHC / Ada / Math_Lib / Polynomials |
1.1 maekawa 1: with unchecked_deallocation;
2: with Generic_Lists;
3: with Graded_Lexicographic_Order; use Graded_Lexicographic_Order;
4:
5: package body Generic_Polynomials is
6:
7: -- REPRESENTATION INVARIANT :
8: -- 1. Only terms with a coefficient different from zero are stored.
9: -- 2. The terms in any polynomial are ordered from high to low degree
10: -- according to the graded lexicographic order.
11:
12: use Ring;
13:
14: -- DATA STRUCTURES :
15:
16: package Term_List is new Generic_Lists(Term);
17: type Poly_Rep is new Term_List.List;
18:
19: procedure free is new unchecked_deallocation(Poly_Rep,Poly);
20:
21: -- AUXILIARY OPERATIONS :
22:
23: procedure Shuffle ( p : in out Poly ) is
24:
25: -- DESCRIPTION :
26: -- Changes the position of the terms in p back to the normal order.
27: -- Needed to guarantee the second representation invariant.
28:
29: res : Poly := Null_Poly;
30:
31: procedure Shuffle_Term ( t : in Term; cont : out boolean ) is
32: begin
33: Add(res,t);
34: cont := true;
35: end Shuffle_Term;
36: procedure Shuffle_Terms is new Visiting_Iterator(Shuffle_Term);
37:
38: begin
39: Shuffle_Terms(p);
40: Clear(p); Copy(res,p); Clear(res);
41: end Shuffle;
42:
43: procedure Append_Copy ( first,last : in out Poly_Rep; t : in Term ) is
44:
45: -- DESCRIPTION :
46: -- Appends a copy of the term to the list.
47:
48: tt : Term;
49:
50: begin
51: Copy(t,tt);
52: Append(first,last,tt);
53: end Append_Copy;
54:
55: -- CONSTRUCTORS :
56:
57: function Create ( n : natural ) return Poly is
58: begin
59: return Create(Ring.Create(n));
60: end Create;
61:
62: function Create ( n : number ) return Poly is
63:
64: t : Term;
65:
66: begin
67: Copy(n,t.cf);
68: return Create(t);
69: end Create;
70:
71: function Create ( t : Term ) return Poly is
72:
73: p : Poly;
74:
75: begin
76: if Equal(t.cf,zero)
77: then p := Null_Poly;
78: else declare
79: tt : Term;
80: begin
81: Copy(t,tt);
82: p := new Poly_Rep;
83: Construct(tt,p.all);
84: end;
85: end if;
86: return p;
87: end Create;
88:
89: procedure Copy ( t1 : in Term; t2 : in out Term ) is
90: begin
91: Clear(t2);
92: Standard_Natural_Vectors.Copy
93: (Standard_Natural_Vectors.Link_to_Vector(t1.dg),
94: Standard_Natural_Vectors.Link_to_Vector(t2.dg));
95: Copy(t1.cf,t2.cf);
96: end Copy;
97:
98: procedure Copy ( p: in Poly_Rep; q : in out Poly_Rep ) is
99:
100: lp,lq : Poly_Rep;
101: t : Term;
102:
103: begin
104: Clear(q);
105: if not Is_Null(p)
106: then lp := p;
107: while not Is_Null(lp) loop
108: t := Head_Of(lp);
109: Append_Copy(q,lq,t);
110: lp := Tail_Of(lp);
111: end loop;
112: end if;
113: end Copy;
114:
115: procedure Copy ( p : in Poly; q : in out Poly ) is
116:
117: l : Poly_Rep;
118:
119: begin
120: Clear(q);
121: if p /= Null_Poly
122: then Copy(p.all,l);
123: q := new Poly_Rep'(l);
124: else q := Null_Poly;
125: end if;
126: end Copy;
127:
128: -- SELECTORS :
129:
130: function Equal ( t1,t2 : Term ) return boolean is
131: begin
132: if not Equal(t1.dg,t2.dg)
133: then return false;
134: else return Equal(t1.cf,t2.cf);
135: end if;
136: end Equal;
137:
138: function Equal ( p,q : Poly_Rep ) return boolean is
139:
140: lp, lq : Poly_Rep;
141:
142: begin
143: lp := p;
144: lq := q;
145: while not Is_Null(lp) and not Is_Null(lq) loop
146: if not Equal(Head_Of(lp),Head_Of(lq))
147: then return false;
148: else lp := Tail_Of(lp);
149: lq := Tail_Of(lq);
150: end if;
151: end loop;
152: if Is_Null(lp) and Is_Null(lq)
153: then return true;
154: else return false;
155: end if;
156: end Equal;
157:
158: function Equal ( p,q : Poly ) return boolean is
159: begin
160: if (p = Null_Poly) and (q = Null_Poly)
161: then return true;
162: elsif (p = Null_Poly) or (q = Null_Poly)
163: then return false;
164: else return Equal(p.all,q.all);
165: end if;
166: end Equal;
167:
168: function Number_Of_Unknowns ( p : Poly ) return natural is
169:
170: temp : Term;
171:
172: begin
173: if (p = Null_Poly) or else Is_Null(p.all)
174: then return 0;
175: else temp := Head_Of(p.all);
176: if temp.dg = null
177: then return 0;
178: else return temp.dg'length;
179: end if;
180: end if;
181: end Number_Of_Unknowns;
182:
183: function Number_Of_Terms ( p : Poly ) return natural is
184: begin
185: if (p = Null_Poly) or else Is_Null(p.all)
186: then return 0;
187: else return Length_Of(p.all);
188: end if;
189: end Number_Of_Terms;
190:
191: function Degree ( p : Poly ) return integer is
192:
193: temp : Term;
194:
195: begin
196: if (p = Null_Poly) or else Is_Null(p.all)
197: then return -1;
198: else temp := Head_Of(p.all);
199: if temp.dg = null
200: then return 0;
201: else return integer(Sum(temp.dg));
202: end if;
203: end if;
204: end Degree;
205:
206: function Degree ( p : Poly; i : integer ) return integer is
207:
208: res : integer := 0;
209:
210: procedure Degree_Term (t : in Term; continue : out boolean) is
211: index,temp : integer;
212: begin
213: if t.dg /= null
214: then index := t.dg'first + i - 1;
215: temp := t.dg(index);
216: if (temp > 0) and then (temp > res)
217: then res := temp;
218: end if;
219: continue := true;
220: end if;
221: end Degree_Term;
222: procedure Degree_Terms is new Visiting_Iterator(process => Degree_Term);
223:
224: begin
225: if p = Null_Poly or else Is_Null(p.all)
226: then return -1;
227: else Degree_Terms(p);
228: return res;
229: end if;
230: end Degree;
231:
232: function "<" ( d1,d2 : Degrees ) return boolean is
233: begin
234: return Standard_Natural_Vectors.Link_to_Vector(d1)
235: < Standard_Natural_Vectors.Link_to_Vector(d2);
236: end "<";
237:
238: function ">" ( d1,d2 : Degrees ) return boolean is
239: begin
240: return Standard_Natural_Vectors.Link_to_Vector(d1)
241: > Standard_Natural_Vectors.Link_to_Vector(d2);
242: end ">";
243:
244: function Coeff ( p : Poly; d : degrees ) return number is
245:
246: l : Poly_Rep;
247: t : term;
248: res : number;
249:
250: begin
251: Copy(zero,res);
252: if p /= Null_Poly
253: then l := p.all;
254: while not Is_Null(l) loop
255: t := Head_Of(l);
256: exit when (t.dg < d);
257: if Equal(t.dg,d)
258: then Copy(t.cf,res);
259: exit;
260: else l := Tail_Of(l);
261: end if;
262: end loop;
263: end if;
264: return res;
265: end Coeff;
266:
267: -- ARITHMETICAL OPERATIONS :
268:
269: function "+" ( p : Poly; t : Term ) return Poly is
270:
271: temp : Poly;
272:
273: begin
274: Copy(p,temp);
275: Add(temp,t);
276: return temp;
277: end "+";
278:
279: function "+" ( t : Term; p : Poly ) return Poly is
280: begin
281: return p+t;
282: end "+";
283:
284: function "+" ( p,q : Poly ) return Poly is
285:
286: temp : Poly;
287:
288: begin
289: Copy(p,temp);
290: Add(temp,q);
291: return temp;
292: end "+";
293:
294: function "+" ( p : Poly ) return Poly is
295:
296: res : Poly;
297:
298: begin
299: Copy(p,res);
300: return res;
301: end "+";
302:
303: function "-" ( p : Poly; t : Term ) return Poly is
304:
305: temp : Poly;
306:
307: begin
308: Copy(p,temp);
309: Sub(temp,t);
310: return temp;
311: end "-";
312:
313: function "-" ( t : Term; p : Poly ) return Poly is
314:
315: temp : Poly;
316:
317: begin
318: temp := Create(t);
319: Sub(temp,p);
320: return temp;
321: end "-";
322:
323: function "-" ( p : Poly ) return Poly is
324:
325: temp : Poly;
326:
327: begin
328: Copy(p,temp);
329: Min(temp);
330: return temp;
331: end "-";
332:
333: function "-" ( p,q : Poly ) return Poly is
334:
335: temp : Poly;
336:
337: begin
338: Copy(p,temp);
339: Sub(temp,q);
340: return temp;
341: end "-";
342:
343: function "*" ( p : Poly; a : number ) return Poly is
344:
345: temp : Poly;
346:
347: begin
348: Copy(p,temp);
349: Mul(temp,a);
350: return temp;
351: end "*";
352:
353: function "*" ( a : number; p : Poly ) return Poly is
354: begin
355: return p*a;
356: end "*";
357:
358: function "*" ( p : Poly; t : Term ) return Poly is
359:
360: temp : Poly;
361:
362: begin
363: Copy(p,temp);
364: Mul(temp,t);
365: return temp;
366: end "*";
367:
368: function "*" ( t : Term; p : Poly ) return Poly is
369: begin
370: return p*t;
371: end "*";
372:
373: function "*" ( p,q : Poly ) return Poly is
374:
375: temp : Poly;
376:
377: begin
378: Copy(p,temp);
379: Mul(temp,q);
380: return temp;
381: end "*";
382:
383: procedure Add ( p : in out Poly; t : in Term ) is
384:
385: l1,l2,temp : Poly_Rep;
386: tt,tp : Term;
387:
388: begin
389: if t.cf /= zero
390: then Copy(t,tt);
391: if p = Null_Poly
392: then p := new Poly_Rep;
393: Construct(tt,p.all);
394: else tp := Head_Of(p.all);
395: if tt.dg > tp.dg
396: then Construct(tt,p.all);
397: elsif Equal(tt.dg,tp.dg)
398: then Add(tp.cf,tt.cf);
399: if tp.cf /= zero
400: then Set_Head(p.all,tp);
401: else Clear(tp);
402: if Is_Null(Tail_Of(p.all))
403: then Term_List.Clear(Term_List.List(p.all));
404: free(p);
405: else Swap_Tail(p.all,l1);
406: Term_List.Clear(Term_List.List(p.all));
407: p.all := l1;
408: end if;
409: end if;
410: Clear(tt);
411: else l1 := p.all;
412: l2 := Tail_Of(l1);
413: while not Is_Null(l2) loop
414: tp := Head_Of(l2);
415: if tt.dg > tp.dg
416: then Construct(tt,temp);
417: Swap_Tail(l1,temp);
418: l1 := Tail_Of(l1);
419: Swap_Tail(l1,temp);
420: return;
421: elsif Equal(tt.dg,tp.dg)
422: then Add(tp.cf,tt.cf);
423: if tp.cf /= zero
424: then Set_Head(l2,tp);
425: else Clear(tp);
426: temp := Tail_Of(l2);
427: Swap_Tail(l1,temp);
428: end if;
429: Clear(tt);
430: return;
431: else l1 := l2;
432: l2 := Tail_Of(l1);
433: end if;
434: end loop;
435: Construct(tt,temp);
436: Swap_Tail(l1,temp);
437: end if;
438: end if;
439: end if;
440: end Add;
441:
442: procedure Add ( p : in out Poly; q : in Poly ) is
443:
444: procedure Add ( t : in Term; continue : out boolean ) is
445: begin
446: Add(p,t);
447: continue := true;
448: end Add;
449: procedure Adds is new Visiting_Iterator(Add);
450:
451: begin
452: Adds(q);
453: end Add;
454:
455: procedure Sub ( p : in out Poly; t : in Term ) is
456:
457: temp : Term;
458:
459: begin
460: Standard_Natural_Vectors.Copy
461: (Standard_Natural_Vectors.Link_to_Vector(t.dg),
462: Standard_Natural_Vectors.Link_to_Vector(temp.dg));
463: temp.cf := -t.cf;
464: Add(p,temp);
465: Standard_Natural_Vectors.Clear
466: (Standard_Natural_Vectors.Link_to_Vector(temp.dg));
467: Clear(temp.cf);
468: end Sub;
469:
470: procedure Sub ( p : in out Poly; q : in Poly ) is
471:
472: temp : Poly := -q;
473:
474: begin
475: Add(p,temp);
476: Clear(temp);
477: end Sub;
478:
479: procedure Min ( p : in out Poly ) is
480:
481: procedure Min ( t : in out Term; continue : out boolean ) is
482: begin
483: Min(t.cf);
484: continue := true;
485: end Min;
486: procedure Min_Terms is new Changing_Iterator (process => Min);
487:
488: begin
489: Min_Terms(p);
490: end Min;
491:
492: procedure Mul ( p : in out Poly; a : in number ) is
493:
494: procedure Mul_Term ( t : in out Term; continue : out boolean ) is
495: begin
496: Mul(t.cf,a);
497: continue := true;
498: end Mul_Term;
499: procedure Mul_Terms is new Changing_Iterator (process => Mul_Term);
500:
501: begin
502: Mul_Terms(p);
503: end Mul;
504:
505: procedure Mul ( p : in out Poly; t : in Term ) is
506:
507: procedure Mul ( tp : in out Term; continue : out boolean ) is
508: begin
509: Standard_Natural_Vectors.Add
510: (Standard_Natural_Vectors.Link_to_Vector(tp.dg),
511: Standard_Natural_Vectors.Link_to_Vector(t.dg));
512: Mul(tp.cf,t.cf);
513: continue := true;
514: end Mul;
515: procedure Mul_Terms is new Changing_Iterator (process => Mul);
516:
517: begin
518: Mul_Terms(p);
519: end Mul;
520:
521: procedure Mul ( p : in out Poly; q : in Poly ) is
522:
523: res : Poly;
524: l : Poly_Rep;
525: t : Term;
526:
527: begin
528: if (q = Null_Poly) or else Is_Null(q.all)
529: then Clear(p);
530: else l := q.all;
531: res := Null_Poly;
532: while not Is_Null(l) loop
533: t := Head_Of(l);
534: declare
535: temp : Poly;
536: begin
537: temp := p * t;
538: Add(res,temp);
539: Clear(temp);
540: end;
541: l := Tail_Of(l);
542: end loop;
543: Copy(res,p); Clear(res);
544: end if;
545: end Mul;
546:
547: procedure Diff ( p : in out Poly; i : in integer ) is
548:
549: procedure Diff_Term ( t : in out Term; continue : out boolean ) is
550:
551: temp : number;
552: index : integer := t.dg'first + i - 1;
553:
554: begin
555: if t.dg(index) = 0
556: then Clear(t);
557: Copy(zero,t.cf);
558: else temp := Create(t.dg(index));
559: Mul(t.cf,temp);
560: Clear(temp);
561: t.dg(index) := t.dg(index) - 1;
562: end if;
563: continue := true;
564: end Diff_Term;
565: procedure Diff_Terms is new Changing_Iterator( process => Diff_Term );
566:
567: begin
568: Diff_Terms(p);
569: end Diff;
570:
571: function Diff ( p : Poly; i : integer ) return Poly is
572:
573: res : Poly;
574:
575: begin
576: Copy(p,res);
577: Diff(res,i);
578: return res;
579: end Diff;
580:
581: -- ITERATORS :
582:
583: procedure Visiting_Iterator ( p : in Poly ) is
584:
585: l : Poly_Rep;
586: temp : Term;
587: continue : boolean;
588:
589: begin
590: if p /= Null_Poly
591: then l := p.all;
592: while not Is_Null(l) loop
593: temp := Head_Of(l);
594: process(temp,continue);
595: exit when not continue;
596: l := Tail_Of(l);
597: end loop;
598: end if;
599: end Visiting_Iterator;
600:
601: procedure Changing_Iterator ( p : in out Poly ) is
602:
603: q,lq,lp : Poly_Rep;
604: t : Term;
605: continue : boolean := true;
606:
607: procedure Copy_append is
608:
609: temp : Term;
610:
611: begin
612: Copy(t,temp);
613: if continue
614: then process(temp,continue);
615: end if;
616: if not Equal(temp.cf,zero)
617: then Append(q,lq,temp);
618: else Clear(temp);
619: end if;
620: Clear(t);
621: end Copy_append;
622:
623: begin
624: if p = Null_Poly
625: then return;
626: else lp := p.all;
627: while not Is_Null(lp) loop
628: t := Head_Of(lp);
629: Copy_append;
630: lp := Tail_Of(lp);
631: end loop;
632: Term_List.Clear(Term_List.List(p.all)); free(p);
633: if Is_Null(q)
634: then p := Null_Poly;
635: else p := new Poly_Rep'(q);
636: end if;
637: Shuffle(p); -- ensure the second representation invariant
638: end if;
639: end Changing_Iterator;
640:
641: -- DESTRUCTORS :
642:
643: procedure Clear ( d : in out Degrees ) is
644: begin
645: Standard_Natural_Vectors.Clear(Standard_Natural_Vectors.Link_to_Vector(d));
646: end Clear;
647:
648: procedure Clear ( t : in out Term ) is
649: begin
650: Clear(t.dg);
651: Clear(t.cf);
652: end Clear;
653:
654: procedure Clear ( p : in out Poly_Rep ) is
655:
656: l : Poly_Rep := p;
657: t : Term;
658:
659: begin
660: if Is_Null(p)
661: then return;
662: else while not Is_Null(l) loop
663: t := Head_Of(l);
664: Clear(t);
665: l := Tail_Of(l);
666: end loop;
667: Term_List.Clear(Term_List.List(p));
668: end if;
669: end Clear;
670:
671: procedure Clear ( p : in out Poly ) is
672: begin
673: if p /= Null_Poly
674: then Clear(p.all);
675: free(p);
676: p := Null_Poly;
677: end if;
678: end Clear;
679:
680: end Generic_Polynomials;