Annotation of OpenXM_contrib/PHC/Ada/Root_Counts/Product/m_homogeneous_bezout_numbers.adb, Revision 1.1
1.1 ! maekawa 1: with Standard_Integer_Vectors; use Standard_Integer_Vectors;
! 2: with Standard_Integer_Matrices; use Standard_Integer_Matrices;
! 3: with Standard_Integer_Linear_Solvers; use Standard_Integer_Linear_Solvers;
! 4: with Standard_Complex_Polynomials; use Standard_Complex_Polynomials;
! 5: with Sets_of_Unknowns; use Sets_of_Unknowns;
! 6: with Degrees_in_Sets_of_Unknowns; use Degrees_in_Sets_of_Unknowns;
! 7:
! 8: package body m_Homogeneous_Bezout_Numbers is
! 9:
! 10: -- UTILITIES :
! 11:
! 12: function Create ( p : Poly_Sys ) return Set is
! 13:
! 14: -- DESCRIPTION :
! 15: -- Returns the set of the unknowns of the polynomial system p.
! 16:
! 17: s : Set := Create(p'length);
! 18:
! 19: begin
! 20: for i in p'range loop
! 21: Add(s,i);
! 22: end loop;
! 23: return s;
! 24: end Create;
! 25:
! 26: function Cardinalities ( z : Partition ) return Vector is
! 27:
! 28: -- DESCRIPTION :
! 29: -- Returns a vector which contains the cardinality of each set.
! 30:
! 31: res : Vector(z'range);
! 32:
! 33: begin
! 34: for i in z'range loop
! 35: res(i) := Extent_Of(z(i));
! 36: end loop;
! 37: return res;
! 38: end Cardinalities;
! 39:
! 40: -- TARGET ROUTINES :
! 41:
! 42: function Total_Degree ( p : Poly_Sys ) return natural is
! 43:
! 44: d : natural := 1;
! 45:
! 46: begin
! 47: for i in p'range loop
! 48: d := d*Degree(p(i));
! 49: end loop;
! 50: return d;
! 51: end Total_Degree;
! 52:
! 53: function Bezout_Number ( p : Poly_Sys; z : Partition ) return natural is
! 54:
! 55: k : constant vector := Cardinalities(z);
! 56: d : constant matrix := Degree_Table(p,z);
! 57:
! 58: begin
! 59: return Per(d,k); -- returns the permanent of the degree table
! 60: end Bezout_Number;
! 61:
! 62: function Bezout_Number ( p : Poly_Sys; z : Partition; max : natural )
! 63: return natural is
! 64:
! 65: -- DESCRIPTION :
! 66: -- Stops when the Bezout number becomes bigger than max.
! 67:
! 68: k : constant vector := Cardinalities(z);
! 69: d : constant matrix := Degree_Table(p,z);
! 70:
! 71: begin
! 72: return Per(d,k,max); -- returns the permanent of the degree table
! 73: end Bezout_Number;
! 74:
! 75: function Bezout_Number ( p : Poly_Sys ) return natural is
! 76:
! 77: s : Set := Create(p);
! 78: res : natural := Total_Degree(p);
! 79:
! 80: procedure Evaluate ( z : in Partition; cont : out boolean ) is
! 81: b : constant natural := Bezout_Number(p,z,res);
! 82: begin
! 83: if b < res
! 84: then res := b;
! 85: end if;
! 86: cont := true;
! 87: end Evaluate;
! 88: procedure Evaluate_Partitions is new Generate_Partitions(Evaluate);
! 89:
! 90: begin
! 91: Evaluate_Partitions(s);
! 92: Clear(s);
! 93: return res;
! 94: end Bezout_Number;
! 95:
! 96: function Bezout_Number ( max : natural; p : Poly_Sys ) return natural is
! 97:
! 98: s : Set := Create(p);
! 99: res : natural := Total_Degree(p);
! 100: cnt : natural := 0;
! 101:
! 102: procedure Evaluate ( z : in Partition; cont : out boolean ) is
! 103: b : constant natural := Bezout_Number(p,z,res);
! 104: begin
! 105: if b < res
! 106: then res := b;
! 107: end if;
! 108: cnt := cnt + 1;
! 109: if cnt < max
! 110: then cont := true;
! 111: else cont := false;
! 112: end if;
! 113: end Evaluate;
! 114: procedure Evaluate_Partitions is new Generate_Partitions(Evaluate);
! 115:
! 116: begin
! 117: Evaluate_Partitions(s);
! 118: Clear(s);
! 119: return res;
! 120: end Bezout_Number;
! 121:
! 122: function Bezout_Number ( p : Poly_Sys; min : natural ) return natural is
! 123:
! 124: s : Set := Create(p);
! 125: res : natural := Total_Degree(p);
! 126:
! 127: procedure Evaluate ( z : in Partition; cont : out boolean ) is
! 128: b : constant natural := Bezout_Number(p,z,res);
! 129: begin
! 130: if b < res
! 131: then res := b;
! 132: end if;
! 133: if res > min
! 134: then cont := true;
! 135: else cont := false;
! 136: end if;
! 137: end Evaluate;
! 138: procedure Evaluate_Partitions is new Generate_Partitions(Evaluate);
! 139:
! 140: begin
! 141: Evaluate_Partitions(s);
! 142: Clear(s);
! 143: return res;
! 144: end Bezout_Number;
! 145:
! 146: function Bezout_Number ( max : natural; p : Poly_Sys; min : natural )
! 147: return natural is
! 148:
! 149: s : Set := Create(p);
! 150: res : natural := Total_Degree(p);
! 151: cnt : natural := 0;
! 152:
! 153: procedure Evaluate ( z : in Partition; cont : out boolean ) is
! 154: b : constant natural := Bezout_Number(p,z,res);
! 155: begin
! 156: if b < res
! 157: then res := b;
! 158: end if;
! 159: cnt := cnt + 1;
! 160: if cnt < max and then res > min
! 161: then cont := true;
! 162: else cont := false;
! 163: end if;
! 164: end Evaluate;
! 165: procedure Evaluate_Partitions is new Generate_Partitions(Evaluate);
! 166:
! 167: begin
! 168: Evaluate_Partitions(s);
! 169: Clear(s);
! 170: return res;
! 171: end Bezout_Number;
! 172:
! 173: procedure Bezout_Number
! 174: ( p : in Poly_Sys; b,m : out natural; z : in out Partition ) is
! 175:
! 176: s : Set := Create(p);
! 177: tdg : constant natural := Total_Degree(p);
! 178: res : natural := tdg;
! 179:
! 180: procedure Evaluate ( nz : in Partition; cont : out boolean ) is
! 181: nb : constant natural := Bezout_Number(p,nz,res);
! 182: begin
! 183: if nb < res
! 184: then res := nb;
! 185: m := nz'length; Clear(z);
! 186: z(nz'range) := Create(nz);
! 187: end if;
! 188: cont := true;
! 189: end Evaluate;
! 190: procedure Evaluate_Partitions is new Generate_Partitions(Evaluate);
! 191:
! 192: begin
! 193: Evaluate_Partitions(s);
! 194: if res = tdg
! 195: then m := 1; z(1) := s;
! 196: else Clear(s);
! 197: end if;
! 198: b := res;
! 199: end Bezout_Number;
! 200:
! 201: procedure Bezout_Number
! 202: ( max : in natural; p : in Poly_Sys; b,m : out natural;
! 203: z : in out Partition ) is
! 204:
! 205: s : Set := Create(p);
! 206: tdg : constant natural := Total_Degree(p);
! 207: res : natural := tdg;
! 208: cnt : natural := 0;
! 209:
! 210: procedure Evaluate ( nz : in Partition; cont : out boolean ) is
! 211: nb : constant natural := Bezout_Number(p,nz,res);
! 212: begin
! 213: if nb < res
! 214: then res := nb;
! 215: m := nz'length; Clear(z);
! 216: z(nz'range) := Create(nz);
! 217: end if;
! 218: cnt := cnt + 1;
! 219: if cnt < max
! 220: then cont := true;
! 221: else cont := false;
! 222: end if;
! 223: end Evaluate;
! 224: procedure Evaluate_Partitions is new Generate_Partitions(Evaluate);
! 225:
! 226: begin
! 227: Evaluate_Partitions(s);
! 228: if res = tdg
! 229: then m := 1; z(1) := s;
! 230: else Clear(s);
! 231: end if;
! 232: b := res;
! 233: end Bezout_Number;
! 234:
! 235: procedure Bezout_Number
! 236: ( p : in Poly_Sys; min : in natural; b,m : out natural;
! 237: z : in out Partition ) is
! 238:
! 239: s : Set := Create(p);
! 240: tdg : constant natural := Total_Degree(p);
! 241: res : natural := tdg;
! 242:
! 243: procedure Evaluate ( nz : in Partition; cont : out boolean ) is
! 244: nb : constant natural := Bezout_Number(p,nz,res);
! 245: begin
! 246: if nb < res
! 247: then res := nb;
! 248: m := nz'length; Clear(z);
! 249: z(nz'range) := Create(nz);
! 250: end if;
! 251: if res > min
! 252: then cont := true;
! 253: else cont := false;
! 254: end if;
! 255: end Evaluate;
! 256: procedure Evaluate_Partitions is new Generate_Partitions(Evaluate);
! 257:
! 258: begin
! 259: Evaluate_Partitions(s);
! 260: if res = tdg
! 261: then m := 1; z(1) := s;
! 262: else Clear(s);
! 263: end if;
! 264: b := res;
! 265: end Bezout_Number;
! 266:
! 267: procedure Bezout_Number
! 268: ( max : in natural; p : in Poly_Sys; min : in natural;
! 269: b,m : out natural; z : in out Partition ) is
! 270:
! 271: s : Set := Create(p);
! 272: tdg : constant natural := Total_Degree(p);
! 273: res : natural := tdg;
! 274: cnt : natural := 0;
! 275:
! 276: procedure Evaluate ( nz : in Partition; cont : out boolean ) is
! 277: nb : constant natural := Bezout_Number(p,nz,res);
! 278: begin
! 279: if nb < res
! 280: then res := nb;
! 281: m := nz'length; Clear(z);
! 282: z(nz'range) := Create(nz);
! 283: end if;
! 284: cnt := cnt + 1;
! 285: if cnt < max and then res > min
! 286: then cont := true;
! 287: else cont := false;
! 288: end if;
! 289: end Evaluate;
! 290: procedure Evaluate_Partitions is new Generate_Partitions(Evaluate);
! 291:
! 292: begin
! 293: Evaluate_Partitions(s);
! 294: if res = tdg
! 295: then m := 1; z(1) := s;
! 296: else Clear(s);
! 297: end if;
! 298: b := res;
! 299: end Bezout_Number;
! 300:
! 301: function Evaluate ( z : partition; m : natural; p : Poly_Sys )
! 302: return natural is
! 303:
! 304: n : natural := p'length;
! 305: d : constant matrix := Degree_Table(p,z);
! 306:
! 307: function Bezout_number
! 308: ( n,m : natural; z: partition; d : matrix ) return natural is
! 309:
! 310: -- DESCRIPTION : the Bezout number is computed
! 311:
! 312: type boolean_array is array ( positive range <> ) of boolean;
! 313: Ii,Iacc : boolean_array(1..n) := (1..n => false);
! 314: b,b_mult : natural;
! 315:
! 316: procedure column ( j,start,number : in natural;
! 317: Ii,Iacc : in out boolean_array );
! 318:
! 319: -- DESCRIPTION : the computation of a term coming from a column
! 320: -- of the degree table;
! 321:
! 322: procedure row ( j : in natural; Ii,Iacc : in out boolean_array );
! 323:
! 324: -- DESCRIPTION : the computation of a row in the degree table
! 325:
! 326: procedure column ( j,start,number : in natural;
! 327: Ii,Iacc : in out boolean_array ) is
! 328: begin
! 329: if number > (n - start + 1)
! 330: then return;
! 331: elsif number = 0
! 332: then row(j,Ii,Iacc);
! 333: else for i in start..n loop
! 334: if not Ii(i)
! 335: then Ii(i) := true; Iacc(i) := true;
! 336: column(j,i+1,number-1,Ii,Iacc);
! 337: Ii(i) := false; Iacc(i) := false;
! 338: end if;
! 339: end loop;
! 340: end if;
! 341: end column;
! 342:
! 343: procedure row ( j : in natural; Ii,Iacc : in out boolean_array ) is
! 344: temp : natural := 1;
! 345: Iacc1 : boolean_array(1..n) := (1..n => false);
! 346: begin
! 347: for k in 1..n loop
! 348: if Iacc(k)
! 349: then temp := temp * d(k,j);
! 350: end if;
! 351: end loop;
! 352: if (j /= m) and (temp /= 0)
! 353: then b_mult := b_mult * temp;
! 354: column(j+1,1,Extent_Of(z(j+1)),Ii,Iacc1);
! 355: b_mult := b_mult / temp;
! 356: elsif j = m
! 357: then temp := temp * b_mult;
! 358: b := b + temp;
! 359: end if;
! 360: end row;
! 361:
! 362: begin
! 363: b := 0; b_mult := 1;
! 364: column(1,1,Extent_Of(z(1)),Ii,Iacc);
! 365: return b;
! 366: end Bezout_number;
! 367:
! 368: begin
! 369: return Bezout_number(n,m,z,d);
! 370: end Evaluate;
! 371:
! 372: procedure PB ( p : in Poly_Sys; b,m : out natural; z : in out Partition ) is
! 373:
! 374: n : natural := p'length;
! 375: wb,b_min : natural;
! 376: wz,z_min : partition(1..n);
! 377: wm,m_min : natural := 0;
! 378:
! 379: procedure pcopy ( p1,p2 : in out partition; p1n,p2n : in out natural ) is
! 380: -- DESCRIPTION : the partition p1 is copied to p2
! 381: begin
! 382: Clear(p2); p2 := Create(p1);
! 383: p2n := p1n;
! 384: end pcopy;
! 385:
! 386: begin
! 387: b_min := Total_Degree(p);
! 388: for i in 1..n loop
! 389: for k in 1..wm loop
! 390: Add(wz(k),i);
! 391: wb := Evaluate(wz,wm,p);
! 392: if (k = 1) or else (wb < b_min)
! 393: then pcopy(wz,z_min,wm,m_min);
! 394: b_min := wb;
! 395: end if;
! 396: Remove(wz(k),i);
! 397: end loop;
! 398: wm := wm + 1;
! 399: wz(wm) := Create(n);
! 400: Add(wz(wm),i);
! 401: wb := Evaluate(wz,wm,p);
! 402: if wb < b_min
! 403: then pcopy(wz,z_min,wm,m_min);
! 404: b_min := wb;
! 405: end if;
! 406: pcopy(z_min,wz,m_min,wm);
! 407: end loop;
! 408: b := b_min;
! 409: m := m_min;
! 410: pcopy(z_min,z,m_min,wm);
! 411: Clear(wz);
! 412: end PB;
! 413:
! 414: end m_Homogeneous_Bezout_Numbers;
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to m_homogeneous_bezout_numbers.adb CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / PHC / Ada / Root_Counts / Product