Annotation of OpenXM_contrib/PHC/Ada/Root_Counts/Stalift/integer_pruning_methods.adb, Revision 1.1
1.1 ! maekawa 1: with Standard_Floating_Numbers; use Standard_Floating_Numbers;
! 2: with Standard_Floating_Matrices;
! 3: with Standard_Integer_Norms; use Standard_Integer_Norms;
! 4: with Standard_Integer_Matrices; use Standard_Integer_Matrices;
! 5: with Standard_Integer_Linear_Solvers; use Standard_Integer_Linear_Solvers;
! 6: with Standard_Integer_Linear_Equalities; use Standard_Integer_Linear_Equalities;
! 7: with Integer_Linear_Inequalities; use Integer_Linear_Inequalities;
! 8: with Floating_Linear_Inequality_Solvers; use Floating_Linear_Inequality_Solvers;
! 9:
! 10: package body Integer_Pruning_Methods is
! 11:
! 12: -- GENERAL AUXILIARIES :
! 13:
! 14: function Convert ( fa : Face_Array ) return Array_of_Lists is
! 15:
! 16: -- DESCRIPTION :
! 17: -- Converts the array of faces into an array of lists by
! 18: -- converting the first element of each list of faces.
! 19:
! 20: res : Array_of_Lists(fa'range);
! 21:
! 22: begin
! 23: for k in fa'range loop
! 24: res(k) := Shallow_Create(fa(k).all);
! 25: end loop;
! 26: return res;
! 27: end Convert;
! 28:
! 29: -- AUXILIARIES FOR THE PRUNING ALGORITHMS :
! 30:
! 31: procedure Initialize ( n : in natural; mat : out Matrix;
! 32: ipvt : out Vector ) is
! 33:
! 34: -- DESCRIPTION :
! 35: -- Initializes the matrix of the equality constraints on the normals
! 36: -- and sets the pivoting vector ipvt to 1..n+1.
! 37:
! 38: begin
! 39: for i in 1..n+1 loop
! 40: ipvt(i) := i;
! 41: for j in 1..n+1 loop
! 42: mat(i,j) := 0;
! 43: end loop;
! 44: end loop;
! 45: end Initialize;
! 46:
! 47: function Number_of_Inequalities
! 48: ( mix : Vector; lifted : Array_of_Lists ) return natural is
! 49:
! 50: -- DESCRIPTION :
! 51: -- Returns the maximal number of inequalities on the inner normal.
! 52:
! 53: res : natural := 0;
! 54:
! 55: begin
! 56: for k in lifted'range loop
! 57: res := res + Length_Of(lifted(k)) - mix(k) - 1;
! 58: end loop;
! 59: return res;
! 60: end Number_of_Inequalities;
! 61:
! 62: procedure Ordered_Inequalities ( k : in natural; mat : in out matrix ) is
! 63:
! 64: -- DESCRIPTION :
! 65: -- Defines k inequalities mat(k,k) - mat(k+1,k) >= 0.
! 66:
! 67: begin
! 68: for i in mat'first(1)..k loop
! 69: for j in mat'range(2) loop
! 70: mat(i,j) := 0;
! 71: end loop;
! 72: mat(i,i) := 1; mat(i,i+1) := -1;
! 73: end loop;
! 74: end Ordered_Inequalities;
! 75:
! 76: procedure Check_and_Update
! 77: ( mic : in Face_Array; lifted : in Array_of_Lists;
! 78: m : in matrix; ipvt : in Vector;
! 79: mixsub,mixsub_last : in out Mixed_Subdivision ) is
! 80:
! 81: -- DESCRIPTION :
! 82: -- Computes the normal to the points in mic, by solving the
! 83: -- linear system defined by m and ipvt. Updates the mixed subdivision
! 84: -- if the computed normal is an inner normal.
! 85:
! 86: v : Vector(m'range(2)) := Solve(m,ipvt);
! 87: pts : Array_of_Lists(mic'range);
! 88:
! 89: begin
! 90: if v(v'last) /= 0
! 91: then Normalize(v);
! 92: if v(v'last) < 0
! 93: then Min(v);
! 94: end if;
! 95: pts := Convert(mic);
! 96: Update(pts,v,mixsub,mixsub_last);
! 97: end if;
! 98: end Check_and_Update;
! 99:
! 100: procedure Create_Equalities
! 101: ( n : in natural; f : in Face; mat,ineq : in matrix;
! 102: ipvt : in Vector; row,rowineq : in natural;
! 103: newmat,newineq : in out matrix; newipvt : in out Vector;
! 104: newrow,newrowineq : in out natural; fail : out boolean ) is
! 105:
! 106: -- DESCRIPTION :
! 107: -- Creates new equalities and uses them to eliminate unknowns in the
! 108: -- matrices of equalities and inequalities. Failure is reported when
! 109: -- a zero row is encountered. On entry, all new* variables must be
! 110: -- initialized with the corresponding *-ones.
! 111:
! 112: shi : vector(1..n+1) := f(f'first).all;
! 113: fl : boolean := false;
! 114: pivot : natural;
! 115:
! 116: begin
! 117: for i in f'first+1..f'last loop
! 118: newrow := newrow + 1;
! 119: for j in f(i)'range loop
! 120: newmat(newrow,j) := f(i)(j) - shi(j);
! 121: end loop;
! 122: Triangulate(newrow,newmat,newipvt,pivot);
! 123: fl := (pivot = 0);
! 124: exit when fl;
! 125: Switch(newrow,pivot,ineq'first,rowineq,newineq);
! 126: --Scale(newrow,newmat);
! 127: end loop;
! 128: fail := fl;
! 129: end Create_Equalities;
! 130:
! 131: function Check_Feasibility
! 132: ( k,m,n : natural; ineq : matrix ) return boolean is
! 133:
! 134: -- DESCRIPTION :
! 135: -- Returns true if -v(n+1) > 0 can be derived, false otherwise.
! 136:
! 137: -- ON ENTRY :
! 138: -- k current unknown that has been eliminated,
! 139: -- for all i <= k : ineq(l,i) = 0, for l in ineq'first..m;
! 140: -- m number of inequalities;
! 141: -- n dimension;
! 142: -- ineq matrix of inequalities.
! 143:
! 144: tableau : matrix(1..n-k+1,1..m+1);
! 145: feasi : boolean;
! 146:
! 147: begin
! 148: if m = 0
! 149: then feasi := false;
! 150: else for i in k+1..n+1 loop
! 151: for j in 1..m loop
! 152: tableau(i-k,j) := ineq(j,i);
! 153: end loop;
! 154: tableau(i-k,m+1) := 0;
! 155: end loop;
! 156: tableau(n-k+1,m+1) := -1;
! 157: Integer_Complementary_Slackness(tableau,feasi);
! 158: end if;
! 159: return feasi;
! 160: end Check_Feasibility;
! 161:
! 162: function New_Check_Feasibility
! 163: ( k,m,n : natural; ineq : matrix; tol : double_float;
! 164: solu : Standard_Floating_Vectors.Vector ) return boolean is
! 165:
! 166: -- DESCRIPTION :
! 167: -- Returns true if -v(n+1) > 0 can be derived, false otherwise.
! 168:
! 169: -- ON ENTRY :
! 170: -- k current unknown that has been eliminated,
! 171: -- for all i <= k : ineq(l,i) = 0, for l in ineq'first..m;
! 172: -- m number of inequalities;
! 173: -- n dimension;
! 174: -- ineq matrix of inequalities.
! 175:
! 176: tableau : Standard_Floating_Matrices.matrix(1..n-k+1,1..m);
! 177: tabsolu : Standard_Floating_Vectors.Vector(1..n-k);
! 178: cols : Vector(tabsolu'range);
! 179: kfail : integer;
! 180: max,tmpabs : double_float;
! 181: -- used for scaling of the last component of the inner normal
! 182:
! 183: begin
! 184: for i in k+1..n loop
! 185: for j in 1..m loop
! 186: tableau(i-k,j) := double_float(ineq(j,i));
! 187: end loop;
! 188: tabsolu(i-k) := solu(i);
! 189: end loop;
! 190: max := 0.0;
! 191: for j in 1..m loop
! 192: tableau(n-k+1,j) := -double_float(ineq(j,n+1));
! 193: tmpabs := abs(tableau(n-k+1,j));
! 194: if tmpabs > max
! 195: then max := tmpabs;
! 196: end if;
! 197: end loop;
! 198: if max > 1.0
! 199: then for j in 1..m loop
! 200: tableau(n-k+1,j) := tableau(n-k+1,j)/max;
! 201: end loop;
! 202: end if;
! 203: Scale(tableau,tol);
! 204: Solve(tableau,tol,tabsolu,kfail,cols);
! 205: return (kfail <= m);
! 206: end New_Check_Feasibility;
! 207:
! 208: procedure Update_Inequalities
! 209: ( k,rowmat1,rowmat2,n : in natural;
! 210: mat : in matrix; ipvt : in vector;
! 211: rowineq : in out natural; ineq : in out matrix;
! 212: lifted : in Array_of_Lists; mic : in out Face_Array ) is
! 213:
! 214: -- DESCRIPTION :
! 215: -- The inequalities will be updated w.r.t. the equality
! 216: -- constraints on the inner normal.
! 217:
! 218: tmp : List;
! 219: pt,shi : Link_to_Vector;
! 220:
! 221: begin
! 222: for i in ineq'first..rowineq loop -- triangulate old inequalities
! 223: for j in rowmat1..rowmat2 loop
! 224: Triangulate(j,mat,i,ineq);
! 225: end loop;
! 226: --Scale(i,ineq);
! 227: end loop;
! 228: tmp := lifted(k); -- create and triangulate
! 229: shi := mic(k)(mic(k)'first); -- new inequalities
! 230: while not Is_Null(tmp) loop
! 231: pt := Head_Of(tmp);
! 232: if not Is_In(mic(k),pt.all)
! 233: then rowineq := rowineq + 1;
! 234: for j in pt'range loop
! 235: ineq(rowineq,j) := pt(j) - shi(j);
! 236: end loop;
! 237: Switch(ipvt,rowineq,ineq);
! 238: for i in 1..rowmat2 loop
! 239: Triangulate(i,mat,rowineq,ineq);
! 240: --Scale(rowineq,ineq);
! 241: end loop;
! 242: end if;
! 243: tmp := Tail_Of(tmp);
! 244: end loop;
! 245: end Update_Inequalities;
! 246:
! 247: procedure Eliminate ( k : in natural; m : in Matrix; tol : in double_float;
! 248: x : in out Standard_Floating_Vectors.Vector ) is
! 249:
! 250: -- DESCRIPTION :
! 251: -- Eliminates the kth component of x by using the matrix.
! 252:
! 253: fac : double_float;
! 254:
! 255: begin
! 256: if abs(x(k)) > tol
! 257: then fac := x(k)/double_float(m(k,k));
! 258: for j in k..x'last loop
! 259: x(j) := x(j) - fac*double_float(m(k,j));
! 260: end loop;
! 261: end if;
! 262: end Eliminate;
! 263:
! 264: procedure Switch ( l,pivot : in integer;
! 265: x : in out Standard_Floating_Vectors.Vector ) is
! 266:
! 267: -- DESCRIPTION :
! 268: -- Applies the pivotation information to the vector x.
! 269:
! 270: tmp : double_float;
! 271:
! 272: begin
! 273: if l /= pivot
! 274: then tmp := x(l); x(l) := x(pivot); x(pivot) := tmp;
! 275: end if;
! 276: end Switch;
! 277:
! 278: procedure Eliminate ( rowmat1,rowmat2 : in natural; mat : in Matrix;
! 279: ipvt : in Vector; tol : in double_float;
! 280: x : in out Standard_Floating_Vectors.Vector ) is
! 281:
! 282: -- DESCRIPTION :
! 283: -- Returns an eliminated solution vector.
! 284:
! 285: begin
! 286: for k in rowmat1..rowmat2 loop
! 287: Switch(k,ipvt(k),x);
! 288: Eliminate(k,mat,tol,x);
! 289: end loop;
! 290: end Eliminate;
! 291:
! 292: -- GENERAL CONSTRUCTOR in several versions :
! 293:
! 294: -- procedure Compute_Mixed_Cells ( ); -- general specification.
! 295:
! 296: -- DESCRIPTION :
! 297: -- Backtrack recursive procedure to enumerate face-face combinations.
! 298:
! 299: -- ON ENTRY :
! 300: -- k index for current point set;
! 301: -- row number of rows already in the matrix;
! 302: -- mat matrix which determines the inner normal;
! 303: -- ipvt contains the pivoting information;
! 304: -- rowineq number of inequality constraints already in ineq;
! 305: -- ineq matrix for the inequality constraints on the
! 306: -- inner normal v, of type <.,v> >= 0;
! 307: -- testsolu test solution to check current set of inequalilities;
! 308: -- mic contains the current selected faces, up to k-1.
! 309:
! 310: -- ON RETURN :
! 311: -- mic updated selected faces;
! 312: -- issupp true if current tuple of faces is supported;
! 313: -- continue indicates whether to continue the creation or not.
! 314:
! 315: -- CONSTRUCTION WITH PRUNING :
! 316:
! 317: procedure Gen1_Create_CS
! 318: ( n : in natural; mix : in Vector; fa : in Array_of_Faces;
! 319: lifted : in Array_of_Lists;
! 320: nbsucc,nbfail : in out Standard_Floating_Vectors.Vector;
! 321: mixsub : out Mixed_Subdivision ) is
! 322:
! 323: res,res_last : Mixed_Subdivision;
! 324: accu : Face_Array(fa'range);
! 325: ma : matrix(1..n+1,1..n+1);
! 326: ipvt : vector(1..n+1);
! 327: ineqrows : natural;
! 328:
! 329: procedure Compute_Mixed_Cells
! 330: ( k,row : in natural; mat : in matrix; ipvt : in vector;
! 331: rowineq : in natural; ineq : in matrix;
! 332: mic : in out Face_Array; continue : out boolean );
! 333:
! 334: -- DESCRIPTION :
! 335: -- Backtrack recursive procedure to enumerate face-face combinations.
! 336:
! 337: procedure Process_Inequalities
! 338: ( k,rowmat1,rowmat2 : in natural;
! 339: mat : in matrix; ipvt : in vector;
! 340: rowineq : in out natural; ineq : in out matrix;
! 341: mic : in out Face_Array; cont : out boolean ) is
! 342:
! 343: -- DESCRIPTION :
! 344: -- Updates inequalities and checks feasibility before proceeding.
! 345:
! 346: begin
! 347: Update_Inequalities(k,rowmat1,rowmat2,n,mat,ipvt,rowineq,ineq,lifted,mic);
! 348: if Check_Feasibility(rowmat2,rowineq,n,ineq)
! 349: then nbfail(k) := nbfail(k) + 1.0;
! 350: cont := true;
! 351: else nbsucc(k) := nbsucc(k) + 1.0;
! 352: Compute_Mixed_Cells(k+1,rowmat2,mat,ipvt,rowineq,ineq,mic,cont);
! 353: end if;
! 354: end Process_Inequalities;
! 355:
! 356: procedure Compute_Mixed_Cells
! 357: ( k,row : in natural; mat : in matrix; ipvt : in vector;
! 358: rowineq : in natural; ineq : in matrix;
! 359: mic : in out Face_Array; continue : out boolean ) is
! 360:
! 361: old : Mixed_Subdivision := res_last;
! 362: cont : boolean := true;
! 363: tmpfa : Faces;
! 364:
! 365: begin
! 366: if k > mic'last
! 367: then Check_and_Update(mic,lifted,mat,ipvt,res,res_last);
! 368: if old /= res_last
! 369: then Process(Head_Of(res_last),continue);
! 370: else continue := true;
! 371: end if;
! 372: else tmpfa := fa(k);
! 373: while not Is_Null(tmpfa) loop -- enumerate faces of kth polytope
! 374: mic(k) := Head_Of(tmpfa);
! 375: declare -- update matrices
! 376: fl : boolean;
! 377: newipvt : vector(ipvt'range) := ipvt;
! 378: newmat : matrix(mat'range(1),mat'range(2)) := mat;
! 379: newineq : matrix(ineq'range(1),ineq'range(2)) := ineq;
! 380: newrow : natural := row;
! 381: newrowineq : natural := rowineq;
! 382: begin
! 383: Create_Equalities(n,mic(k),mat,ineq,ipvt,row,rowineq,
! 384: newmat,newineq,newipvt,newrow,newrowineq,fl);
! 385: if fl
! 386: then nbfail(k) := nbfail(k) + 1.0;
! 387: else Process_Inequalities(k,row+1,newrow,newmat,newipvt,
! 388: newrowineq,newineq,mic,cont);
! 389: end if;
! 390: end;
! 391: exit when not cont;
! 392: tmpfa := Tail_Of(tmpfa);
! 393: end loop;
! 394: continue := cont;
! 395: end if;
! 396: end Compute_Mixed_Cells;
! 397:
! 398: begin
! 399: Initialize(n,ma,ipvt);
! 400: ineqrows := Number_of_Inequalities(mix,lifted);
! 401: declare
! 402: ineq : matrix(1..ineqrows,1..n+1);
! 403: cont : boolean;
! 404: begin
! 405: ineq(1,1) := 0;
! 406: Compute_Mixed_Cells(accu'first,0,ma,ipvt,0,ineq,accu,cont);
! 407: end;
! 408: mixsub := res;
! 409: end Gen1_Create_CS;
! 410:
! 411: procedure Create_CS
! 412: ( n : in natural; mix : in Vector; fa : in Array_of_Faces;
! 413: lifted : in Array_of_Lists;
! 414: nbsucc,nbfail : in out Standard_Floating_Vectors.Vector;
! 415: mixsub : out Mixed_Subdivision ) is
! 416:
! 417: res,res_last : Mixed_Subdivision;
! 418: accu : Face_Array(fa'range);
! 419: ma : matrix(1..n+1,1..n+1);
! 420: ipvt : vector(1..n+1);
! 421: ineqrows : natural;
! 422:
! 423: procedure Compute_Mixed_Cells
! 424: ( k,row : in natural; mat : in Matrix; ipvt : in Vector;
! 425: rowineq : in natural; ineq : in Matrix;
! 426: mic : in out Face_Array );
! 427:
! 428: -- DESCRIPTION :
! 429: -- Backtrack recursive procedure to enumerate face-face combinations.
! 430:
! 431: procedure Process_Inequalities
! 432: ( k,rowmat1,rowmat2 : in natural;
! 433: mat : in matrix; ipvt : in vector;
! 434: rowineq : in out natural; ineq : in out matrix;
! 435: mic : in out Face_Array ) is
! 436:
! 437: -- DESCRIPTION :
! 438: -- Updates inequalities and checks feasibility before proceeding.
! 439:
! 440: begin
! 441: Update_Inequalities(k,rowmat1,rowmat2,n,mat,ipvt,rowineq,ineq,lifted,mic);
! 442: if Check_Feasibility(rowmat2,rowineq,n,ineq)
! 443: then nbfail(k) := nbfail(k) + 1.0;
! 444: else nbsucc(k) := nbsucc(k) + 1.0;
! 445: Compute_Mixed_Cells(k+1,rowmat2,mat,ipvt,rowineq,ineq,mic);
! 446: end if;
! 447: end Process_Inequalities;
! 448:
! 449: procedure Compute_Mixed_Cells
! 450: ( k,row : in natural; mat : in Matrix; ipvt : in Vector;
! 451: rowineq : in natural; ineq : in Matrix;
! 452: mic : in out Face_Array ) is
! 453:
! 454: tmpfa : Faces;
! 455:
! 456: begin
! 457: if k > mic'last
! 458: then Check_and_Update(mic,lifted,mat,ipvt,res,res_last);
! 459: else tmpfa := fa(k);
! 460: while not Is_Null(tmpfa) loop -- enumerate faces of kth polytype
! 461: mic(k) := Head_Of(tmpfa);
! 462: declare -- update matrices
! 463: fl : boolean;
! 464: newipvt : vector(ipvt'range) := ipvt;
! 465: newmat : matrix(mat'range(1),mat'range(2)) := mat;
! 466: newineq : matrix(ineq'range(1),ineq'range(2)) := ineq;
! 467: newrow : natural := row;
! 468: newrowineq : natural := rowineq;
! 469: begin
! 470: Create_Equalities(n,mic(k),mat,ineq,ipvt,row,rowineq,newmat,
! 471: newineq,newipvt,newrow,newrowineq,fl);
! 472: if fl
! 473: then nbfail(k) := nbfail(k) + 1.0;
! 474: else Process_Inequalities
! 475: (k,row+1,newrow,newmat,newipvt,newrowineq,newineq,mic);
! 476: end if;
! 477: end;
! 478: tmpfa := Tail_Of(tmpfa);
! 479: end loop;
! 480: end if;
! 481: end Compute_Mixed_Cells;
! 482:
! 483: begin
! 484: Initialize(n,ma,ipvt);
! 485: ineqrows := Number_of_Inequalities(mix,lifted);
! 486: declare
! 487: ineq : matrix(1..ineqrows,1..n+1);
! 488: begin
! 489: ineq(1,1) := 0;
! 490: Compute_Mixed_Cells(accu'first,0,ma,ipvt,0,ineq,accu);
! 491: end;
! 492: mixsub := res;
! 493: end Create_CS;
! 494:
! 495: procedure New_Create_CS
! 496: ( n : in natural; mix : in Vector; fa : in Array_of_Faces;
! 497: lifted : in Array_of_Lists;
! 498: nbsucc,nbfail : in out Standard_Floating_Vectors.Vector;
! 499: mixsub : out Mixed_Subdivision ) is
! 500:
! 501: res,res_last : Mixed_Subdivision;
! 502: accu : Face_Array(fa'range);
! 503: ma : matrix(1..n+1,1..n+1);
! 504: ipvt : vector(1..n+1);
! 505: tol : constant double_float := double_float(n)*10.0**(-11);
! 506: ineqrows : natural;
! 507:
! 508: procedure Compute_Mixed_Cells
! 509: ( k,row : in natural; mat : in matrix; ipvt : in vector;
! 510: rowineq : in natural; ineq : in matrix;
! 511: testsolu : in Standard_Floating_Vectors.Vector;
! 512: mic : in out Face_Array );
! 513:
! 514: -- DESCRIPTION :
! 515: -- Backtrack recursive procedure to enumerate face-face combinations.
! 516:
! 517: procedure Process_Inequalities
! 518: ( k,rowmat1,rowmat2 : in natural;
! 519: mat : in matrix; ipvt : in vector;
! 520: rowineq : in out natural; ineq : in out matrix;
! 521: testsolu : in Standard_Floating_Vectors.Vector;
! 522: mic : in out Face_Array ) is
! 523:
! 524: -- DESCRIPTION :
! 525: -- Updates the inequalities and checks feasibility before proceeding.
! 526:
! 527: tmp : List;
! 528: pt,shi : Link_to_Vector;
! 529: newsolu : Standard_Floating_Vectors.Vector(testsolu'range) := testsolu;
! 530:
! 531: begin
! 532: Update_Inequalities(k,rowmat1,rowmat2,n,mat,ipvt,rowineq,ineq,
! 533: lifted,mic);
! 534: Eliminate(rowmat1,rowmat2,mat,ipvt,tol,newsolu);
! 535: if New_Check_Feasibility(rowmat2,rowineq,n,ineq,tol,newsolu)
! 536: then nbfail(k) := nbfail(k) + 1.0;
! 537: else nbsucc(k) := nbsucc(k) + 1.0;
! 538: Compute_Mixed_Cells(k+1,rowmat2,mat,ipvt,rowineq,ineq,testsolu,mic);
! 539: end if;
! 540: end Process_Inequalities;
! 541:
! 542: procedure Compute_Mixed_Cells
! 543: ( k,row : in natural; mat : in matrix; ipvt : in vector;
! 544: rowineq : in natural; ineq : in matrix;
! 545: testsolu : in Standard_Floating_Vectors.Vector;
! 546: mic : in out Face_Array ) is
! 547:
! 548: -- DESCRIPTION :
! 549: -- Backtrack recursive procedure to enumerate face-face combinations.
! 550:
! 551: tmpfa : Faces;
! 552:
! 553: begin
! 554: if k > mic'last
! 555: then Check_and_Update(mic,lifted,mat,ipvt,res,res_last);
! 556: else tmpfa := fa(k);
! 557: while not Is_Null(tmpfa) loop -- enumerate faces of kth polytope
! 558: mic(k) := Head_Of(tmpfa);
! 559: declare -- update matrices
! 560: fl : boolean;
! 561: newipvt : vector(ipvt'range) := ipvt;
! 562: newmat : matrix(mat'range(1),mat'range(2)) := mat;
! 563: newineq : matrix(ineq'range(1),ineq'range(2)) := ineq;
! 564: newrow : natural := row;
! 565: newrowineq : natural := rowineq;
! 566: begin
! 567: Create_Equalities(n,mic(k),mat,ineq,ipvt,row,rowineq,newmat,
! 568: newineq,newipvt,newrow,newrowineq,fl);
! 569: if fl
! 570: then nbfail(k) := nbfail(k) + 1.0;
! 571: else Process_Inequalities(k,row+1,newrow,newmat,newipvt,
! 572: newrowineq,newineq,testsolu,mic);
! 573: end if;
! 574: end;
! 575: tmpfa := Tail_Of(tmpfa);
! 576: end loop;
! 577: end if;
! 578: end Compute_Mixed_Cells;
! 579:
! 580: begin
! 581: Initialize(n,ma,ipvt);
! 582: ineqrows := Number_of_Inequalities(mix,lifted);
! 583: declare
! 584: ineq : matrix(1..ineqrows,1..n+1);
! 585: solu : Standard_Floating_Vectors.Vector(1..n+1);
! 586: begin
! 587: ineq(1,1) := 0;
! 588: solu := (solu'range => 0.0);
! 589: Compute_Mixed_Cells(accu'first,0,ma,ipvt,0,ineq,solu,accu);
! 590: end;
! 591: mixsub := res;
! 592: end New_Create_CS;
! 593:
! 594: -- AUXILIARIES FOR THE CONSTRAINT PRUNING :
! 595:
! 596: function Is_Supported ( f : Face; normal : Vector; supp : integer )
! 597: return boolean is
! 598:
! 599: ip : integer;
! 600:
! 601: begin
! 602: for i in f'range loop
! 603: ip := f(i).all*normal;
! 604: if ip /= supp
! 605: then return false;
! 606: end if;
! 607: end loop;
! 608: return true;
! 609: end Is_Supported;
! 610:
! 611: function Is_Supported ( f : Face; k : natural; normals,suppvals : List )
! 612: return boolean is
! 613:
! 614: -- DESCRIPTION :
! 615: -- Returns true if there exists a normal for which the inner product
! 616: -- with each point in the face equals the corresponding supporting value
! 617: -- of the kth component.
! 618:
! 619: tmpnor : List := normals;
! 620: tmpsup : List := suppvals;
! 621: support : integer;
! 622:
! 623: begin
! 624: while not Is_Null(tmpnor) loop
! 625: support := Head_Of(tmpsup)(k);
! 626: if Is_Supported(f,Head_Of(tmpnor).all,support)
! 627: then return true;
! 628: else tmpnor := Tail_Of(tmpnor);
! 629: tmpsup := Tail_Of(tmpsup);
! 630: end if;
! 631: end loop;
! 632: return false;
! 633: end Is_Supported;
! 634:
! 635: procedure Update ( mic : in Mixed_Cell; normals,suppvals : in out List ) is
! 636:
! 637: -- DESCRIPTION :
! 638: -- Updates the list of normals and supporting values with the information
! 639: -- of a mixed cell.
! 640:
! 641: normal : Link_to_Vector := new Vector'(mic.nor.all);
! 642: suppval : Link_to_Vector := new Vector(mic.pts'range);
! 643:
! 644: begin
! 645: Construct(normal,normals);
! 646: for i in suppval'range loop
! 647: suppval(i) := normal*Head_Of(mic.pts(i));
! 648: end loop;
! 649: Construct(suppval,suppvals);
! 650: end Update;
! 651:
! 652: procedure Create_CCS
! 653: ( n : in natural; mix : in Vector; fa : in Array_of_Faces;
! 654: lifted : in Array_of_Lists;
! 655: nbsucc,nbfail : in out Standard_Floating_Vectors.Vector;
! 656: normals,suppvals : in out List;
! 657: mixsub : out Mixed_Subdivision ) is
! 658:
! 659: res,res_last : Mixed_Subdivision;
! 660: accu : Face_Array(fa'range);
! 661: ma : matrix(1..n+1,1..n+1);
! 662: ipvt : vector(1..n+1);
! 663: ineqrows : natural;
! 664:
! 665: procedure Compute_Mixed_Cells
! 666: ( k,row : in natural; mat : in matrix; ipvt : in vector;
! 667: rowineq : in natural; ineq : in matrix;
! 668: mic : in out Face_Array; issupp : in out boolean );
! 669:
! 670: -- DESCRIPTION :
! 671: -- Backtrack recursive procedure to enumerate face-face combinations.
! 672:
! 673: procedure Process_Inequalities
! 674: ( k,rowmat1,rowmat2 : in natural;
! 675: mat : in matrix; ipvt : in vector;
! 676: rowineq : in out natural; ineq : in out matrix;
! 677: mic : in out Face_Array; issupp : in out boolean ) is
! 678:
! 679: -- DESCRIPTION :
! 680: -- Updates inequalities and checks feasibility before proceeding.
! 681:
! 682: begin
! 683: Update_Inequalities(k,rowmat1,rowmat2,n,mat,ipvt,rowineq,ineq,lifted,mic);
! 684: if issupp
! 685: then issupp := Is_Supported(mic(k),k,normals,suppvals);
! 686: end if;
! 687: if not issupp and then Check_Feasibility(rowmat2,rowineq,n,ineq)
! 688: then nbfail(k) := nbfail(k) + 1.0;
! 689: else nbsucc(k) := nbsucc(k) + 1.0;
! 690: Compute_Mixed_Cells(k+1,rowmat2,mat,ipvt,rowineq,ineq,mic,issupp);
! 691: end if;
! 692: end Process_Inequalities;
! 693:
! 694: procedure Compute_Mixed_Cells
! 695: ( k,row : in natural; mat : in matrix; ipvt : in vector;
! 696: rowineq : in natural; ineq : in matrix;
! 697: mic : in out Face_Array; issupp : in out boolean ) is
! 698:
! 699: -- DESCRIPTION :
! 700: -- Backtrack recursive procedure to enumerate face-face combinations.
! 701:
! 702: old : Mixed_Subdivision;
! 703: tmpfa : Faces;
! 704:
! 705: begin
! 706: if k > mic'last
! 707: then old := res_last;
! 708: Check_and_Update(mic,lifted,mat,ipvt,res,res_last);
! 709: if old /= res_last
! 710: then Update(Head_Of(res_last),normals,suppvals);
! 711: end if;
! 712: else tmpfa := fa(k);
! 713: while not Is_Null(tmpfa) loop -- enumerate faces of kth polytope
! 714: mic(k) := Head_Of(tmpfa);
! 715: declare -- update matrices
! 716: fl : boolean;
! 717: newipvt : vector(ipvt'range) := ipvt;
! 718: newmat : matrix(mat'range(1),mat'range(2)) := mat;
! 719: newineq : matrix(ineq'range(1),ineq'range(2)) := ineq;
! 720: newrow : natural := row;
! 721: newrowineq : natural := rowineq;
! 722: begin
! 723: Create_Equalities(n,mic(k),mat,ineq,ipvt,row,rowineq,newmat,
! 724: newineq,newipvt,newrow,newrowineq,fl);
! 725: if fl
! 726: then nbfail(k) := nbfail(k) + 1.0;
! 727: else Process_Inequalities(k,row+1,newrow,newmat,newipvt,
! 728: newrowineq,newineq,mic,issupp);
! 729: end if;
! 730: end;
! 731: tmpfa := Tail_Of(tmpfa);
! 732: end loop;
! 733: end if;
! 734: end Compute_Mixed_Cells;
! 735:
! 736: begin
! 737: Initialize(n,ma,ipvt);
! 738: ineqrows := Number_of_Inequalities(mix,lifted);
! 739: declare
! 740: ineq : matrix(1..ineqrows,1..n+1);
! 741: supported : boolean := true;
! 742: begin
! 743: ineq(1,1) := 0;
! 744: Compute_Mixed_Cells(accu'first,0,ma,ipvt,0,ineq,accu,supported);
! 745: end;
! 746: mixsub := res;
! 747: end Create_CCS;
! 748:
! 749: end Integer_Pruning_Methods;
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to integer_pruning_methods.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 / Stalift