Annotation of OpenXM_contrib/PHC/Ada/Math_Lib/Supports/ts_fvector.adb, Revision 1.1.1.1
1.1 maekawa 1: with text_io,integer_io; use text_io,integer_io;
2: with Standard_Floating_Numbers; use Standard_Floating_Numbers;
3: with Standard_Random_Numbers; use Standard_Random_Numbers;
4: with Standard_Integer_Vectors;
5: with Standard_Floating_Vectors;
6: with Standard_Integer_VecVecs;
7: with Standard_Integer_VecVecs_io; use Standard_Integer_VecVecs_io;
8: with Standard_Floating_VecVecs;
9: with Standard_Floating_VecVecs_io; use Standard_Floating_VecVecs_io;
10: with Standard_Random_VecVecs; use Standard_Random_VecVecs;
11: with Face_Cardinalities; use Face_Cardinalities;
12:
13: procedure ts_fvector is
14:
15: -- DESCRIPTION :
16: -- Computes the f-vector of a polytope and checks the Euler-Poincare formula.
17:
18: procedure Write ( f : in Standard_Integer_Vectors.Vector ) is
19:
20: -- DESCRIPTION :
21: -- Writes the f-vector on screen.
22:
23: n : constant natural := f'last;
24:
25: begin
26: put_line("The f-vector : ");
27: put(" f(-1) : "); put(f(-1)); new_line;
28: put(" #vertices : "); put(f(0)); new_line;
29: put(" #edges : "); put(f(1)); new_line;
30: for i in 2..(n-1) loop
31: put(" #"); put(i,1); put("-faces : "); put(f(i)); new_line;
32: end loop;
33: put(" f("); put(n,1); put(") : "); put(f(n)); new_line;
34: end Write;
35:
36: function Euler_Poincare ( f : Standard_Integer_Vectors.Vector )
37: return integer is
38:
39: -- DESCRIPTION :
40: -- Computes the alternating sum: sum_{i in f'range} (-1)^(i) f(i).
41:
42: sum : integer := 0;
43: pos : boolean := false;
44:
45: begin
46: for i in f'range loop
47: if pos
48: then sum := sum + f(i);
49: else sum := sum - f(i);
50: end if;
51: pos := not pos;
52: end loop;
53: return sum;
54: end Euler_Poincare;
55:
56: procedure Integer_Interactive_Testing is
57:
58: n,m : natural;
59: tol : constant double_float := 10.0**(-12);
60: sum : integer;
61:
62: begin
63: put("Give the dimension n : "); get(n);
64: put("Give the number of points that span the polytope : "); get(m);
65: declare
66: pts : Standard_Integer_VecVecs.VecVec(1..m);
67: f : Standard_Integer_Vectors.Vector(-1..n);
68: begin
69: put("Give "); put(m,1); put(" "); put(n,1);
70: put_line("-dimensional integer vectors :");
71: get(n,pts);
72: put_line("Counting vertices, edges , .., facets, ...");
73: f := fvector(pts); sum := Euler_Poincare(f);
74: Write(f);
75: put("The result of the Euler-Poincare formula : "); put(sum,1);
76: if sum /= 0
77: then put_line(" BUG DISCOVERED !!!");
78: else put_line(" OK.");
79: end if;
80: end;
81: end Integer_Interactive_Testing;
82:
83: procedure Floating_Interactive_Testing is
84:
85: n,m : natural;
86: tol : constant double_float := 10.0**(-12);
87: sum : integer;
88:
89: begin
90: put("Give the dimension n : "); get(n);
91: put("Give the number of points that span the polytope : "); get(m);
92: declare
93: pts : Standard_Floating_VecVecs.VecVec(1..m);
94: f : Standard_Integer_Vectors.Vector(-1..n);
95: begin
96: put("Give "); put(m,1); put(" "); put(n,1);
97: put_line("-dimensional floating point vectors :");
98: get(n,pts);
99: put_line("Counting vertices, edges , .., facets, ...");
100: f := fvector(pts); sum := Euler_Poincare(f);
101: Write(f);
102: put("The result of the Euler-Poincare formula : "); put(sum,1);
103: if sum /= 0
104: then put_line(" BUG DISCOVERED !!!");
105: else put_line(" OK.");
106: end if;
107: end;
108: end Floating_Interactive_Testing;
109:
110: function Floating_Random_Polytope
111: ( n,m : natural ) return Standard_Floating_VecVecs.VecVec is
112:
113: -- DESCRIPTION :
114: -- Returns m randomly chosen n-dimensional floating-point vectors.
115:
116: res : Standard_Floating_VecVecs.VecVec(1..m) := Random_VecVec(n,m);
117:
118: begin
119: return res;
120: end Floating_Random_Polytope;
121:
122: function Integer_Random_Polytope
123: ( n,m : natural; lower,upper : integer )
124: return Standard_Integer_VecVecs.VecVec is
125:
126: -- DESCRIPTION :
127: -- Returns m randomly chosen n-dimensional vectors,
128: -- with integer entries between lower and upper.
129:
130: res : Standard_Integer_VecVecs.VecVec(1..m);
131: done : boolean;
132:
133: function Is_In ( i : integer ) return boolean is
134:
135: -- DESCRIPTION :
136: -- Returns true if the ith vector already occurs in res(1..i-1).
137:
138: use Standard_Integer_Vectors;
139:
140: begin
141: for j in 1..i-1 loop
142: if Equal(res(j).all,res(i).all)
143: then return true;
144: end if;
145: end loop;
146: return false;
147: end Is_In;
148:
149: begin
150: for i in 1..m loop
151: res(i) := new Standard_Integer_Vectors.Vector(1..n);
152: done := false;
153: while not done loop
154: for j in 1..n loop
155: res(i)(j) := Random(lower,upper);
156: end loop;
157: done := not Is_In(i);
158: end loop;
159: end loop;
160: return res;
161: end Integer_Random_Polytope;
162:
163: procedure Floating_Automatic_Testing is
164:
165: n,m,times,cnt : natural;
166: tol : constant double_float := 10.0**(-12);
167: sum : integer;
168: bug : boolean := false;
169:
170: begin
171: put("Give the dimension n : "); get(n);
172: put("Give the number of points that span the polytope : "); get(m);
173: put("Give the number of testing cycles : "); get(times);
174: declare
175: pts : Standard_Floating_VecVecs.VecVec(1..m);
176: f : Standard_Integer_Vectors.Vector(-1..n);
177: begin
178: for i in 1..times loop
179: cnt := i;
180: pts := Floating_Random_Polytope(n,m);
181: f := fvector(pts); sum := Euler_Poincare(f);
182: Write(f);
183: put("The result of the Euler-Poincare formula : "); put(sum,1);
184: if sum /= 0
185: then put_line(" BUG DISCOVERED !!!"); bug := true;
186: put_line("The generated random configuration is");
187: put(pts);
188: else put_line(" OK."); bug := false;
189: end if;
190: Standard_Floating_VecVecs.Clear(pts);
191: exit when bug;
192: end loop;
193: end;
194: if not bug
195: then put("No bugs found, with "); put(times,1);
196: put_line(" generated cases tested.");
197: put("Dimension : "); put(n,1);
198: put(" and #points : "); put(m,1); put_line(".");
199: else put("Bug found at case "); put(cnt,1); put_line(".");
200: end if;
201: end Floating_Automatic_Testing;
202:
203: procedure Integer_Automatic_Testing is
204:
205: n,m,times,cnt : natural;
206: tol : constant double_float := 10.0**(-12);
207: sum : integer;
208: bug : boolean := false;
209: lower,upper : integer;
210:
211: begin
212: put("Give the dimension n : "); get(n);
213: put("Give the number of points that span the polytope : "); get(m);
214: put("Give lower bound on the entries : "); get(lower);
215: put("Give upper bound on the entries : "); get(upper);
216: put("Give the number of testing cycles : "); get(times);
217: declare
218: pts : Standard_Integer_VecVecs.VecVec(1..m);
219: f : Standard_Integer_Vectors.Vector(-1..n);
220: begin
221: for i in 1..times loop
222: cnt := i;
223: pts := Integer_Random_Polytope(n,m,lower,upper);
224: f := fvector(pts); sum := Euler_Poincare(f);
225: Write(f);
226: put("The result of the Euler-Poincare formula : "); put(sum,1);
227: if sum /= 0
228: then put_line(" BUG DISCOVERED !!!"); bug := true;
229: put_line("The generated random configuration is");
230: put(pts);
231: else put_line(" OK."); bug := false;
232: end if;
233: Standard_Integer_VecVecs.Clear(pts);
234: exit when bug;
235: end loop;
236: end;
237: if not bug
238: then put("No bugs found, with "); put(times,1);
239: put_line(" generated cases tested.");
240: put("Dimension : "); put(n,1);
241: put(" and #points : "); put(m,1); put_line(".");
242: else put("Bug found at case "); put(cnt,1); put_line(".");
243: end if;
244: end Integer_Automatic_Testing;
245:
246: procedure Main is
247:
248: ans : character;
249:
250: begin
251: new_line;
252: put_line("Testing the face enumerators by computing f-vectors.");
253: loop
254: new_line;
255: put_line("Choose one of the following : ");
256: put_line(" 0. Exit this program. ");
257: put_line(" 1. f-vector of given integer polytope. ");
258: put_line(" 2. f-vector of given floating polytope. ");
259: put_line(" 3. f-vector of random integer polytope. ");
260: put_line(" 4. f-vector of random floating polytope. ");
261: put("Type 0,1,2,3, or 4 to select : "); get(ans);
262: exit when (ans = '0');
263: case ans is
264: when '1' => Integer_Interactive_Testing;
265: when '2' => Floating_Interactive_Testing;
266: when '3' => Integer_Automatic_Testing;
267: when '4' => Floating_Automatic_Testing;
268: when others => null;
269: end case;
270: end loop;
271: end Main;
272:
273: begin
274: Main;
275: end ts_fvector;
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>