Annotation of OpenXM_contrib/PHC/Ada/Schubert/symbolic_minor_equations.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_Natural_Vectors;
! 4: with Matrix_Indeterminates;
! 5: with Bracket_Polynomials; use Bracket_Polynomials;
! 6: with Straightening_Syzygies; use Straightening_Syzygies;
! 7:
! 8: package body Symbolic_Minor_Equations is
! 9:
! 10: -- AUXILIARIES TO Minor_Equations :
! 11:
! 12: function Substitute ( b,minor : Bracket ) return Bracket is
! 13:
! 14: -- DESCRIPTION :
! 15: -- The ith entry in the bracket b is replaced by minor(b(i)).
! 16:
! 17: res : Bracket(b'range);
! 18: start : integer;
! 19:
! 20: begin
! 21: if b(b'first) = 0
! 22: then start := b'first+1; -- bracket is a coefficient bracket
! 23: res(res'first) := 0; -- result is also a coefficient bracket
! 24: else start := b'first; -- second bracket in expansion
! 25: end if;
! 26: for i in start..b'last loop
! 27: res(i) := minor(b(i));
! 28: end loop;
! 29: return res;
! 30: end Substitute;
! 31:
! 32: function Substitute ( bm : Bracket_Monomial; minor : Bracket )
! 33: return Bracket_Monomial is
! 34:
! 35: -- DESCRIPTION :
! 36: -- Substitutes the entries in the brackets according to the minor.
! 37:
! 38: -- REQUIRED : bm is a quadratic monomial.
! 39:
! 40: res : Bracket_Monomial;
! 41: lb : Link_to_Bracket;
! 42: first : boolean := true;
! 43:
! 44: procedure Substitute_Bracket ( b : in Bracket; continue : out boolean ) is
! 45:
! 46: sb : Bracket(b'range) := Substitute(b,minor);
! 47:
! 48: begin
! 49: if first
! 50: then lb := new Bracket'(sb); -- save to preserve order
! 51: first := false;
! 52: else res := sb*lb.all;
! 53: end if;
! 54: continue := true;
! 55: end Substitute_Bracket;
! 56: procedure Substitute_Brackets is
! 57: new Enumerate_Brackets(Substitute_Bracket);
! 58:
! 59: begin
! 60: Substitute_Brackets(bm);
! 61: Clear(lb);
! 62: return res;
! 63: end Substitute;
! 64:
! 65: function Substitute ( bt : Bracket_Term; minor : Bracket )
! 66: return Bracket_Term is
! 67:
! 68: -- DESCRIPTION :
! 69: -- The ith entry in every bracket of the term according to the minor.
! 70:
! 71: res : Bracket_Term;
! 72:
! 73: begin
! 74: res.coeff := bt.coeff;
! 75: res.monom := Substitute(bt.monom,minor);
! 76: return res;
! 77: end Substitute;
! 78:
! 79: function Substitute ( bp : Bracket_Polynomial; minor : Bracket )
! 80: return Bracket_Polynomial is
! 81:
! 82: -- DESCRIPTION :
! 83: -- The labels in the brackets of bp are replaced by those in the minor.
! 84: -- The bracket polynomial represents the Laplace expansion of that minor.
! 85:
! 86: res : Bracket_Polynomial;
! 87:
! 88: procedure Substitute_Term ( t : in Bracket_Term; continue : out boolean ) is
! 89:
! 90: st : Bracket_Term := Substitute(t,minor);
! 91:
! 92: begin
! 93: Add(res,st);
! 94: Clear(st);
! 95: continue := true;
! 96: end Substitute_Term;
! 97: procedure Substitute_Terms is new Enumerate_Terms(Substitute_Term);
! 98:
! 99: begin
! 100: Substitute_Terms(bp);
! 101: return res;
! 102: end Substitute;
! 103:
! 104: -- AUXILIARIES TO Expanded_Minor :
! 105:
! 106: function Subtract ( b : Bracket; i : natural ) return Bracket is
! 107:
! 108: -- DESCRIPTION :
! 109: -- Returns a smaller bracket with the ith entry removed.
! 110:
! 111: res : Bracket(b'first..b'last-1);
! 112:
! 113: begin
! 114: res(b'first..(i-1)) := b(b'first..(i-1));
! 115: res(i..res'last) := b((i+1)..b'last);
! 116: return res;
! 117: end Subtract;
! 118:
! 119: function General_Expanded_Minor ( m : Standard_Complex_Poly_Matrices.Matrix;
! 120: b : Bracket ) return Poly is
! 121:
! 122: -- DESCRIPTION :
! 123: -- This function treats the case for Expanded_Minor when b'length > 2.
! 124:
! 125: res : Poly := Null_Poly;
! 126: sig : integer;
! 127: submin,acc : Poly;
! 128:
! 129: begin
! 130: if b'last mod 2 = 0
! 131: then sig := -1;
! 132: else sig := +1;
! 133: end if;
! 134: for i in b'range loop
! 135: if m(b(i),b'last) /= Null_Poly
! 136: then submin := Expanded_Minor(m,Subtract(b,i));
! 137: if submin /= Null_Poly
! 138: then acc := m(b(i),b'last)*submin;
! 139: if sig > 0
! 140: then Add(res,acc);
! 141: else Sub(res,acc);
! 142: end if;
! 143: Clear(acc);
! 144: end if;
! 145: Clear(submin);
! 146: end if;
! 147: sig := -sig;
! 148: end loop;
! 149: return res;
! 150: end General_Expanded_Minor;
! 151:
! 152: procedure Purify ( p : in out Poly ) is
! 153:
! 154: -- DESCRIPTION :
! 155: -- Eliminates terms that have coefficients less that 10.0**(-10).
! 156:
! 157: tol : constant double_float := 10.0**(-10);
! 158: res : Poly := Null_Poly;
! 159:
! 160: procedure Scan_Term ( t : in Term; cont : out boolean ) is
! 161: begin
! 162: if AbsVal(t.cf) > tol
! 163: then Add(res,t);
! 164: end if;
! 165: cont := true;
! 166: end Scan_Term;
! 167: procedure Scan_Terms is new Visiting_Iterator(Scan_Term);
! 168:
! 169: begin
! 170: Scan_Terms(p);
! 171: Clear(p);
! 172: p := res;
! 173: end Purify;
! 174:
! 175: -- TARGET ROUTINES :
! 176:
! 177: function Schubert_Pattern ( n : natural; b1,b2 : Bracket )
! 178: return Standard_Complex_Poly_Matrices.Matrix is
! 179:
! 180: res : Standard_Complex_Poly_Matrices.Matrix(1..n,b1'range);
! 181: d : constant natural := b1'last;
! 182:
! 183: begin
! 184: for i in 1..n loop
! 185: for j in b1'range loop
! 186: if ((i < b2(j)) or (i > n+1-b1(b1'last+1-j)))
! 187: then res(i,j) := Null_Poly;
! 188: else res(i,j) := Matrix_Indeterminates.Monomial(n,d,i,j);
! 189: end if;
! 190: end loop;
! 191: end loop;
! 192: return res;
! 193: end Schubert_Pattern;
! 194:
! 195: function Localization_Pattern ( n : natural; top,bottom : Bracket )
! 196: return Standard_Complex_Poly_Matrices.Matrix is
! 197:
! 198: p : constant natural := top'length;
! 199: res : Standard_Complex_Poly_Matrices.Matrix(1..n,1..p);
! 200:
! 201: begin
! 202: for i in res'range(1) loop
! 203: for j in res'range(2) loop
! 204: if i < top(j) or i > bottom(j)
! 205: then res(i,j) := Null_Poly;
! 206: else res(i,j) := Matrix_Indeterminates.Monomial(n,p,i,j);
! 207: end if;
! 208: end loop;
! 209: end loop;
! 210: return res;
! 211: end Localization_Pattern;
! 212:
! 213: -- SYMBOLIC REPRESENTATION OF THE EQUATIONS :
! 214:
! 215: function Maximal_Minors ( n,m : natural ) return Bracket_Monomial is
! 216:
! 217: res : Bracket_Monomial;
! 218: first : boolean := true;
! 219: accu : Bracket(1..m);
! 220:
! 221: procedure Enumerate_Brackets ( k,start : natural ) is
! 222:
! 223: -- DESCRIPTION :
! 224: -- Enumerate all m-brackets with entries between 1 and n,
! 225: -- beginning at start, and currently filling in the kth entry.
! 226:
! 227: begin
! 228: if k > m
! 229: then if first
! 230: then res := Create(accu); first := false;
! 231: else Multiply(res,accu);
! 232: end if;
! 233: else for i in start..n loop
! 234: accu(k) := i;
! 235: Enumerate_Brackets(k+1,i+1);
! 236: end loop;
! 237: end if;
! 238: end Enumerate_Brackets;
! 239:
! 240: begin
! 241: Enumerate_Brackets(1,1);
! 242: return res;
! 243: end Maximal_Minors;
! 244:
! 245: function Minor_Equations
! 246: ( m,d : natural; bm : Bracket_Monomial ) return Bracket_System is
! 247:
! 248: res : Bracket_System(0..Number_of_Brackets(bm));
! 249: gen : Bracket_Polynomial := Laplace_Expansion(m,d);
! 250: cnt : natural := 0;
! 251:
! 252: procedure Substitute_Bracket ( b : in Bracket; continue : out boolean ) is
! 253: begin
! 254: cnt := cnt + 1;
! 255: res(cnt) := Substitute(gen,b);
! 256: continue := true;
! 257: end Substitute_Bracket;
! 258: procedure Substitute_Brackets is
! 259: new Enumerate_Brackets(Substitute_Bracket);
! 260:
! 261: begin
! 262: res(0) := gen;
! 263: Substitute_Brackets(bm);
! 264: return res;
! 265: end Minor_Equations;
! 266:
! 267: function Expanded_Minor ( m : Standard_Complex_Poly_Matrices.Matrix;
! 268: b : Bracket ) return Poly is
! 269:
! 270: res : Poly := Null_Poly;
! 271: acc : Poly;
! 272:
! 273: begin
! 274: if b'length = 1
! 275: then Copy(m(b(1),1),res);
! 276: elsif b'length = 2
! 277: then if (m(b(1),1) /= Null_Poly) and (m(b(2),2) /= Null_Poly)
! 278: then res := m(b(1),1)*m(b(2),2);
! 279: end if;
! 280: if (m(b(2),1) /= Null_Poly) and (m(b(1),2) /= Null_Poly)
! 281: then acc := m(b(2),1)*m(b(1),2);
! 282: Sub(res,acc);
! 283: Clear(acc);
! 284: end if;
! 285: else res := General_Expanded_Minor(m,b);
! 286: end if;
! 287: Purify(res);
! 288: return res;
! 289: end Expanded_Minor;
! 290:
! 291: function Extend_Zero_Lifting ( p : Poly ) return Poly is
! 292:
! 293: res : Poly := Null_Poly;
! 294:
! 295: procedure Extend_Term ( t : in Term; continue : out boolean ) is
! 296:
! 297: et : Term;
! 298:
! 299: begin
! 300: et.dg := new Standard_Natural_Vectors.Vector(t.dg'first..t.dg'last+1);
! 301: et.dg(t.dg'range) := t.dg.all;
! 302: et.dg(et.dg'last) := 0;
! 303: et.cf := t.cf;
! 304: Add(res,et);
! 305: Clear(et.dg);
! 306: continue := true;
! 307: end Extend_Term;
! 308: procedure Extend_Terms is new Visiting_Iterator(Extend_Term);
! 309:
! 310: begin
! 311: Extend_Terms(p);
! 312: return res;
! 313: end Extend_Zero_Lifting;
! 314:
! 315: end Symbolic_Minor_Equations;
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to symbolic_minor_equations.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