Annotation of OpenXM_contrib/PHC/Ada/Root_Counts/Product/multi_homogeneous_start_systems.adb, Revision 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: with Standard_Complex_Vectors;
! 6: with Standard_Complex_Matrices; use Standard_Complex_Matrices;
! 7: with Standard_Complex_Linear_Solvers; use Standard_Complex_Linear_Solvers;
! 8: with Standard_Complex_Polynomials; use Standard_Complex_Polynomials;
! 9:
! 10: with Sets_of_Unknowns; use Sets_of_Unknowns;
! 11: with Degrees_in_Sets_of_Unknowns; use Degrees_in_Sets_of_Unknowns;
! 12: with Partitions_of_Sets_of_Unknowns; use Partitions_of_Sets_of_Unknowns;
! 13: with Degree_Structure; use Degree_Structure;
! 14: with Random_Product_System;
! 15:
! 16: package body Multi_Homogeneous_Start_Systems is
! 17:
! 18: procedure Add_Hyperplanes ( i : in natural; p : in Poly;
! 19: z : in partition; m : in natural;
! 20: dg : in Standard_Natural_Vectors.Vector) is
! 21: -- DESCRIPTION :
! 22: -- This routine adds hyperplanes to the Random Product System
! 23: -- according to the partition and the degrees of the polynomial
! 24:
! 25: -- ON ENTRY :
! 26: -- i the number of the polynomial in the system;
! 27: -- p the i-the polynomial of the system;
! 28: -- z a partition of the set of unknowns of p;
! 29: -- m the number of sets in z;
! 30: -- dg the degrees of the polynomial in the sets of z,
! 31: -- for j in 1..m : dg(j) = Degree(p,z(j)).
! 32:
! 33: n : natural := Number_Of_Unknowns(p);
! 34: h : Standard_Complex_Vectors.Vector(0..n);
! 35:
! 36: begin
! 37: for k in 1..m loop
! 38: for l in 1..dg(k) loop
! 39: h(0) := random1;
! 40: for j in 1..n loop
! 41: if Is_In(z(k),j)
! 42: then h(j) := Random1;
! 43: else h(j) := Create(0.0);
! 44: end if;
! 45: end loop;
! 46: Random_Product_System.Add_Hyperplane(i,h);
! 47: end loop;
! 48: end loop;
! 49: end Add_Hyperplanes;
! 50:
! 51: procedure Construct_Random_Product_System ( p : in Poly_Sys ) is
! 52:
! 53: -- DESCRIPTION :
! 54: -- This procedure constructs a random product system by
! 55: -- finding a good partition for each equation of the system p.
! 56:
! 57: n : natural := p'length;
! 58: m : natural;
! 59: begin
! 60: if Degree_Structure.Empty
! 61: then Degree_Structure.Create(p);
! 62: end if;
! 63: for i in p'range loop
! 64: m := Degree_Structure.Get(i);
! 65: declare
! 66: z : Partition(1..m);
! 67: dg : Standard_Natural_Vectors.Vector(1..m);
! 68: begin
! 69: Degree_Structure.Get(i,z,dg);
! 70: Add_Hyperplanes(i,p(i),z,m,dg);
! 71: Clear(z);
! 72: end;
! 73: end loop;
! 74: end Construct_Random_Product_System;
! 75:
! 76: procedure RPQ ( p : in Poly_Sys; q : out Poly_Sys;
! 77: sols : in out Solution_List; nl : out natural) is
! 78:
! 79: n : natural := p'length;
! 80:
! 81: begin
! 82: Random_Product_System.Init(n);
! 83: Construct_Random_Product_System(p);
! 84: q := Random_Product_System.Polynomial_System;
! 85: Random_Product_System.Solve(sols,nl);
! 86: Degree_Structure.Clear;
! 87: Random_Product_System.Clear;
! 88: end RPQ;
! 89:
! 90: procedure Solve ( i,n : in natural; sols : in out Solution_List;
! 91: a : in out Matrix;
! 92: b : in out Standard_Complex_Vectors.Vector;
! 93: acc : in out partition ) is
! 94: begin
! 95: if i > n
! 96: then declare
! 97: aa : Matrix(a'range(1),a'range(2));
! 98: bb : Standard_Complex_Vectors.Vector(b'range);
! 99: rcond : double_float;
! 100: ipvt : Standard_Natural_Vectors.Vector(b'range);
! 101: begin
! 102: for k in a'range(1) loop
! 103: for l in a'range(2) loop
! 104: aa(k,l) := a(k,l);
! 105: end loop;
! 106: bb(k) := b(k);
! 107: end loop;
! 108: lufco(aa,n,ipvt,rcond);
! 109: -- put("rcond : "); put(rcond); new_line;
! 110: if rcond + 1.0 /= 1.0
! 111: then lusolve(aa,n,ipvt,bb);
! 112: declare
! 113: s : Solution(n);
! 114: begin
! 115: s.t := Create(0.0);
! 116: s.m := 1;
! 117: s.v := bb;
! 118: Add(sols,s);
! 119: end;
! 120: end if;
! 121: exception
! 122: when numeric_error => return;
! 123: end;
! 124: else declare
! 125: h : Standard_Complex_Vectors.Vector(0..n);
! 126: count : natural := 0;
! 127: z : Partition(1..n);
! 128: m : natural;
! 129: d : Standard_Natural_Vectors.Vector(1..n);
! 130: begin
! 131: Degree_Structure.Get(i,z,d);
! 132: m := Degree_Structure.Get(i);
! 133: for j1 in 1..m loop
! 134: if Degree_Structure.Admissible(acc,i-1,z(j1))
! 135: then acc(i) := Create(z(j1));
! 136: for j2 in 1..d(j1) loop
! 137: count := count + 1;
! 138: h := Random_Product_System.Get_Hyperplane(i,count);
! 139: b(i) := -h(0);
! 140: for k in 1..n loop
! 141: a(i,k) := h(k);
! 142: end loop;
! 143: Solve(i+1,n,sols,a,b,acc);
! 144: end loop;
! 145: Clear(acc(i));
! 146: else count := count + d(j1);
! 147: end if;
! 148: end loop;
! 149: Clear(z);
! 150: end;
! 151: end if;
! 152: end Solve;
! 153:
! 154: procedure Solve_Start_System
! 155: ( n : in natural; sols : in out Solution_List ) is
! 156:
! 157: m : Matrix(1..n,1..n);
! 158: v : Standard_Complex_Vectors.Vector(1..n);
! 159: acc : Partition(1..n);
! 160:
! 161: begin
! 162: for i in 1..n loop
! 163: for j in 1..n loop
! 164: m(i,j) := Create(0.0);
! 165: end loop;
! 166: v(i) := Create(0.0);
! 167: end loop;
! 168: Solve(1,n,sols,m,v,acc);
! 169: end Solve_Start_System;
! 170:
! 171: procedure GBQ ( p : in Poly_Sys; q : out Poly_Sys;
! 172: sols : in out Solution_List ) is
! 173:
! 174: n : natural := p'length;
! 175:
! 176: begin
! 177: Random_Product_System.Init(n);
! 178: Construct_Random_Product_System(p);
! 179: q := Random_Product_System.Polynomial_System;
! 180: Solve_Start_System(n,sols);
! 181: Degree_Structure.Clear;
! 182: Random_Product_System.Clear;
! 183: end GBQ;
! 184:
! 185: end Multi_Homogeneous_Start_Systems;
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to multi_homogeneous_start_systems.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