Annotation of OpenXM_contrib/PHC/Ada/Homotopy/homogenization.adb, Revision 1.1.1.1
1.1 maekawa 1: with Standard_Floating_Numbers; use Standard_Floating_Numbers;
2: with Standard_Complex_Numbers; use Standard_Complex_Numbers;
3: with Standard_Random_Numbers; use Standard_Random_Numbers;
4: with Standard_Natural_Vectors; use Standard_Natural_Vectors;
5:
6: package body Homogenization is
7:
8: function Homogeneous_Part ( p : Poly ) return Poly is
9:
10: res : Poly := Null_Poly;
11: d : integer := Degree(p);
12:
13: procedure Homogeneous_Term ( t : in Term; continue : out boolean ) is
14: begin
15: if Sum(t.dg) = d
16: then Add(res,t);
17: continue := true;
18: else continue := false;
19: end if;
20: end Homogeneous_Term;
21: procedure Homogeneous_Terms is new Visiting_Iterator(Homogeneous_Term);
22:
23: begin
24: Homogeneous_Terms(p);
25: return res;
26: end Homogeneous_Part;
27:
28: function Homogeneous_Part ( p : Poly_Sys ) return Poly_Sys is
29:
30: res : Poly_Sys(p'range);
31:
32: begin
33: for i in p'range loop
34: res(i) := Homogeneous_Part(p(i));
35: end loop;
36: return res;
37: end Homogeneous_Part;
38:
39: function Real_Random_Hyperplane ( n : natural ) return Poly is
40:
41: res : Poly;
42: t : Term;
43:
44: begin
45: t.dg := new Standard_Natural_Vectors.Vector'(1..n => 0);
46: t.cf := Create(Random);
47: res := Create(t);
48: Standard_Natural_Vectors.Clear
49: (Standard_Natural_Vectors.Link_to_Vector(t.dg));
50: for i in 1..n loop
51: t.dg := new Standard_Natural_Vectors.Vector'(1..n => 0);
52: t.dg(i) := 1;
53: t.cf := Create(Random);
54: Add(res,t);
55: Standard_Natural_Vectors.Clear
56: (Standard_Natural_Vectors.Link_to_Vector(t.dg));
57: end loop;
58: return res;
59: end Real_Random_Hyperplane;
60:
61: function Complex_Random_Hyperplane ( n : natural ) return Poly is
62:
63: res : Poly;
64: t : Term;
65:
66: begin
67: t.dg := new Standard_Natural_Vectors.Vector'(1..n => 0);
68: t.cf := Random;
69: res := Create(t);
70: Standard_Natural_Vectors.Clear
71: (Standard_Natural_Vectors.Link_to_Vector(t.dg));
72: for i in 1..n loop
73: t.dg := new Standard_Natural_Vectors.Vector'(1..n => 0);
74: t.dg(i) := 1;
75: t.cf := Random;
76: Add(res,t);
77: Standard_Natural_Vectors.Clear
78: (Standard_Natural_Vectors.Link_to_Vector(t.dg));
79: end loop;
80: return res;
81: end Complex_Random_Hyperplane;
82:
83: function Standard_Hyperplane ( n,i : natural ) return Poly is
84:
85: -- DESCRIPTION : Returns x_i - 1.
86:
87: res : Poly;
88: t : Term;
89:
90: begin
91: t.dg := new Standard_Natural_Vectors.Vector'(1..n => 0);
92: t.cf := Create(-1.0);
93: res := Create(t);
94: Standard_Natural_Vectors.Clear
95: (Standard_Natural_Vectors.Link_to_Vector(t.dg));
96: t.dg := new Standard_Natural_Vectors.Vector'(1..n => 0);
97: t.dg(i) := 1;
98: t.cf := Create(1.0);
99: Add(res,t);
100: Standard_Natural_Vectors.Clear
101: (Standard_Natural_Vectors.Link_to_Vector(t.dg));
102: return res;
103: end Standard_Hyperplane;
104:
105: procedure Construct_Real_Random_Hyperplanes
106: ( s : in out Poly_Sys; m : natural ) is
107:
108: -- DESCRIPTION :
109: -- the polynomial system s will be filled with polynomials in m unknowns,
110: -- with real coefficients.
111:
112: begin
113: for i in s'range loop
114: Clear(s(i));
115: s(i) := Real_Random_Hyperplane(m);
116: end loop;
117: end Construct_Real_Random_Hyperplanes;
118:
119: procedure Construct_Complex_Random_Hyperplanes
120: ( s : in out Poly_Sys; m : natural ) is
121:
122: -- DESCRIPTION :
123: -- The polynomial system s will be filled with polynomials in m unknowns,
124: -- with complex coefficients.
125:
126: begin
127: for i in s'range loop
128: Clear(s(i));
129: s(i) := Complex_Random_Hyperplane(m);
130: end loop;
131: end Construct_Complex_Random_Hyperplanes;
132:
133: procedure Construct_Standard_Hyperplanes ( s : in out Poly_Sys ) is
134:
135: -- DESCRIPTION :
136: -- for i in s'range : s(i) := x_i - 1.
137:
138: n : natural := s'length;
139:
140: begin
141: for i in s'range loop
142: Clear(s(i));
143: s(i) := Standard_Hyperplane(n,i);
144: end loop;
145: end Construct_Standard_Hyperplanes;
146:
147: procedure Enlarge_Before ( p : in out Poly; m : in natural ) is
148:
149: -- DESCRIPTION :
150: -- To each term t of p, m additional zero entries will be inserted to t.dg
151:
152: procedure Enlarge_Term ( t : in out Term; continue : out boolean ) is
153:
154: d : Degrees := new Standard_Natural_Vectors.Vector(1..(t.dg'last+m));
155:
156: begin
157: for i in 1..m loop
158: d(i) := 0;
159: end loop;
160: for i in t.dg'range loop
161: d(i+m) := t.dg(i);
162: end loop;
163: Standard_Natural_Vectors.Clear
164: (Standard_Natural_Vectors.Link_to_Vector(t.dg));
165: t.dg := d;
166: continue := true;
167: end Enlarge_Term;
168: procedure Enlarge_Terms is new Changing_Iterator(Enlarge_Term);
169:
170: begin
171: Enlarge_Terms(p);
172: end Enlarge_Before;
173:
174: procedure Enlarge_After ( p : in out Poly; m : in natural ) is
175:
176: -- DESCRIPTION :
177: -- To each term t of p, m additional zero entries will be added to t.dg
178:
179: procedure Enlarge_Term ( t : in out Term; continue : out boolean ) is
180:
181: d : Degrees := new Standard_Natural_Vectors.Vector(1..(t.dg'last+m));
182:
183: begin
184: for i in t.dg'range loop
185: d(i) := t.dg(i);
186: end loop;
187: for i in (t.dg'last+1)..m loop
188: d(i) := 0;
189: end loop;
190: Standard_Natural_Vectors.Clear
191: (Standard_Natural_Vectors.Link_to_Vector(t.dg));
192: t.dg := d;
193: continue := true;
194: end Enlarge_Term;
195: procedure Enlarge_Terms is new Changing_Iterator(Enlarge_Term);
196:
197: begin
198: Enlarge_Terms(p);
199: end Enlarge_After;
200:
201: function Add_Equations ( s1 : Poly_Sys; s2 : Poly_Sys ) return Poly_Sys is
202:
203: n1 : natural := s1'last - s1'first + 1;
204: n2 : natural := s2'last - s2'first + 1;
205: res : Poly_Sys(1..(n1+n2));
206: m : natural;
207:
208: begin
209: for i in 1..n1 loop
210: Copy(s1(i),res(i));
211: m := Number_Of_Unknowns(res(i));
212: if m < (n1+n2)
213: then m := n1 + n2 - m;
214: Enlarge_After(res(i),m);
215: end if;
216: end loop;
217: for i in 1..n2 loop
218: Copy(s2(i),res(n1+i));
219: m := Number_Of_Unknowns(res(n1+i));
220: if m < (n1+n2)
221: then m := n1 + n2 - m;
222: Enlarge_Before(res(n1+i),m);
223: end if;
224: end loop;
225: return res;
226: end Add_Equations;
227:
228: function Add_Equation ( s : Poly_Sys; p : Poly ) return Poly_Sys is
229:
230: n : natural := s'last - s'first + 1;
231: m : natural;
232: res : Poly_Sys(1..(n+1));
233:
234: begin
235: for i in 1..n loop
236: Copy(s(i),res(i));
237: m := Number_Of_Unknowns(res(i));
238: if m < n+1
239: then m := n + 1 - m;
240: Enlarge_After(res(i),m);
241: end if;
242: end loop;
243: m := Number_Of_Unknowns(p);
244: if m < (n+1)
245: then m := n + 1 - m;
246: Enlarge_Before(res(n+1),m);
247: end if;
248: return res;
249: end Add_Equation;
250:
251: function Add_Random_Hyperplanes
252: ( s : Poly_Sys; m : natural; re : boolean ) return Poly_Sys is
253:
254: s2 : Poly_Sys(1..m);
255: n : natural := s'last - s'first + 1;
256: res : Poly_Sys(1..(m+n));
257:
258: begin
259: if re
260: then Construct_Real_Random_Hyperplanes(s2,m+n);
261: else Construct_Complex_Random_Hyperplanes(s2,m+n);
262: end if;
263: res := Add_Equations(s,s2);
264: Clear(s2);
265: return res;
266: end Add_Random_Hyperplanes;
267:
268: function Add_Standard_Hyperplanes
269: ( s : Poly_Sys; m : natural ) return Poly_Sys is
270:
271: n : natural := s'length;
272: res : Poly_Sys(1..(n+m));
273: s2 : Poly_Sys(1..m);
274:
275: begin
276: Construct_Standard_Hyperplanes(s2);
277: res := Add_Equations(s,s2);
278: Clear(s2);
279: return res;
280: end Add_Standard_Hyperplanes;
281:
282: end Homogenization;
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to homogenization.adb CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / PHC / Ada / Homotopy