Annotation of OpenXM_contrib/PHC/Ada/Schubert/ts_detrock.adb, Revision 1.1.1.1
1.1 maekawa 1: with text_io,integer_io; use text_io,integer_io;
2: with Characters_and_Numbers; use Characters_and_Numbers;
3: with Communications_with_User; use Communications_with_User;
4: with Standard_Floating_Numbers; use Standard_Floating_Numbers;
5: with Standard_Random_Numbers; use Standard_Random_Numbers;
6: with Standard_Floating_Matrices; use Standard_Floating_Matrices;
7: with Standard_Complex_Poly_Systems; use Standard_Complex_Poly_Systems;
8: with Standard_Complex_Polynomials; use Standard_Complex_Polynomials;
9: with Osculating_Planes; use Osculating_Planes;
10: with Standard_Complex_Numbers; use Standard_Complex_Numbers;
11: with Standard_Complex_Poly_SysFun; use Standard_Complex_Poly_SysFun;
12: with Standard_Complex_Poly_Systems_io; use Standard_Complex_Poly_Systems_io;
13: with Total_Degree_Start_Systems; use Total_Degree_Start_Systems;
14: with Sets_of_Unknowns; use Sets_of_Unknowns;
15: with Partitions_of_Sets_Of_Unknowns; use Partitions_of_Sets_of_Unknowns;
16: with Partitions_of_Sets_Of_Unknowns_io; use Partitions_of_Sets_of_Unknowns_io;
17: with m_Homogeneous_Bezout_Numbers; use m_Homogeneous_Bezout_Numbers;
18: with m_Homogeneous_Start_Systems; use m_Homogeneous_Start_Systems;
19: with Lists_of_Integer_Vectors; use Lists_of_Integer_Vectors;
20: with Lists_of_Integer_Vectors_io; use Lists_of_Integer_Vectors_io;
21: with Power_Lists; use Power_Lists;
22: with Triangulations,Triangulations_io; use Triangulations,Triangulations_io;
23: with Dynamic_Triangulations; use Dynamic_Triangulations;
24: with Matrix_Indeterminates;
25: with Bracket_Expansions; use Bracket_Expansions;
26: with SAGBI_Homotopies; use SAGBI_Homotopies;
27:
28: procedure ts_detrock is
29:
30: -- DESCRIPTION :
31: -- Generates (m,p)-system and performs root counting.
32:
33: m,p : natural;
34: ans : character;
35: file : file_type;
36:
37: function Determinant_System ( m,p : natural ) return Poly_Sys is
38:
39: res : Poly_Sys(1..m*p);
40: lp : Poly;
41: s : double_float := Random;
42: inc : constant double_float := 2.0/double_float(m*p);
43: mat : Matrix(1..m+p,1..m);
44:
45: begin
46: Matrix_Indeterminates.Initialize_Symbols(m+p,p);
47: lp := Lifted_Localized_Laplace_Expansion(m+p,p);
48: for i in res'range loop
49: s := s+inc;
50: mat := Orthogonal_Basis(m+p,m,s);
51: res(i) := Intersection_Condition(mat,lp);
52: if s > 2.0
53: then s := s - 2.0;
54: end if;
55: end loop;
56: return res;
57: end Determinant_System;
58:
59: procedure Count_Roots ( file : in file_type; h : in Poly_Sys;
60: m,p : in natural; title_banner : in string ) is
61:
62: function Minimum ( a,b : natural ) return natural is
63: begin
64: if a <= b
65: then return a;
66: else return b;
67: end if;
68: end Minimum;
69:
70: function Construct_Partition ( m,p : natural ) return Partition is
71:
72: min_mp : constant natural := Minimum(m,p);
73: z : Partition(1..min_mp);
74: cnt : natural := 0;
75:
76: begin
77: for i in z'range loop
78: z(i) := Create(m*p);
79: end loop;
80: if m <= p
81: then for i in 1..m loop
82: for j in 1..p loop
83: cnt := cnt+1;
84: Add(z(i),cnt);
85: end loop;
86: end loop;
87: else cnt := 1;
88: for i in 1..p loop
89: for j in 1..m loop
90: Add(z(i),cnt);
91: cnt := cnt+p;
92: if cnt > m*p
93: then cnt := cnt-m*p;
94: end if;
95: end loop;
96: end loop;
97: end if;
98: return z;
99: end Construct_Partition;
100:
101: procedure Multi_Homogeneous_Bound ( f : in Poly_Sys ) is
102:
103: b,nz : natural;
104: -- z : Partition(p'range);
105: min_mp : constant natural := Minimum(m,p);
106: z : Partition(1..min_mp) := Construct_Partition(m,p);
107:
108: begin
109: -- PB(f,b,nz,z);
110: nz := z'last;
111: b := Bezout_Number(f,z);
112: put(file,nz,1); put(file,"-homogeneous Bezout number : ");
113: put(file,b,1); new_line(file);
114: put(file," with partition "); put(file,z); new_line(file);
115: Clear(z);
116: end Multi_Homogeneous_Bound;
117:
118: procedure Apply_Root_Counts ( f : in Poly_Sys; cmpvol : in boolean ) is
119:
120: d : natural := Total_Degree(f);
121: sup,lifted,lifted_last : List;
122: t : Triangulation;
123: vol : natural;
124:
125: begin
126: new_line(file);
127: put_line(file,"ROOT COUNTS : ");
128: new_line(file);
129: put(file,"total degree : "); put(file,d,1); new_line(file);
130: Multi_Homogeneous_Bound(f);
131: if cmpvol
132: then sup := Create(f(f'first));
133: Dynamic_Lifting(sup,false,false,0,lifted,lifted_last,t);
134: vol := Volume(t);
135: put(file,"mixed volume : "); put(file,vol,1); new_line(file);
136: new_line(file);
137: put_line(file,"The lifted support : ");
138: new_line(file);
139: put(file,lifted);
140: -- new_line(file);
141: -- put_line(file,"The regular triangulation : ");
142: -- new_line(file);
143: -- put(file,p'length,t,vol);
144: Clear(t); Clear(sup); Clear(lifted);
145: end if;
146: end Apply_Root_Counts;
147:
148: procedure Main is
149:
150: target,start : Poly_Sys(h'range);
151:
152: begin
153: target := Eval(h,Create(1.0),m*p+1);
154: put(file,target'length,target);
155: new_line(file);
156: put_line(file,title_banner);
157: Apply_Root_Counts(target,false);
158: start := Eval(h,Create(0.0),m*p+1);
159: new_line(file);
160: put(file,start'length,start);
161: new_line(file);
162: put_line(file,"TITLE : start system in SAGBI homotopy.");
163: Apply_Root_Counts(start,true);
164: Clear(target); Clear(start);
165: end Main;
166:
167: begin
168: Main;
169: end Count_Roots;
170:
171: begin
172: new_line;
173: put_line("Performing root counts on determinantal (m,p)-systems.");
174: loop
175: new_line;
176: put_line("Reading the name of the output file.");
177: Read_Name_and_Create_File(file);
178: put("Give m : "); get(m);
179: put("Give p : "); get(p);
180: declare
181: title : constant string := "TITLE : all " & Convert(p)
182: & "-planes that intersect " & Convert(m*p)
183: & " random real osculating " & Convert(m)
184: & "-planes.";
185: begin
186: -- put(file,"Determinantal ("); put(file,m,1); put(file,",");
187: -- put(file,p,1); put_line(file,")-system :");
188: Count_Roots(file,Determinant_System(m,p),m,p,title);
189: end;
190: Close(file);
191: put("Do you want other systems to test ? (y/n) ");
192: Ask_Yes_or_No(ans);
193: exit when (ans /= 'y');
194: end loop;
195: end ts_detrock;
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to ts_detrock.adb CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / PHC / Ada / Schubert