Annotation of OpenXM_contrib/PHC/Ada/Root_Counts/Symmetry/orbits_of_solutions.adb, Revision 1.1
1.1 ! maekawa 1: with unchecked_deallocation;
! 2: with Standard_Complex_Numbers; use Standard_Complex_Numbers;
! 3: with Standard_Complex_Vectors;
! 4: with Standard_Complex_Norms_Equals; use Standard_Complex_Norms_Equals;
! 5: with Permute_Operations; use Permute_Operations;
! 6:
! 7: package body Orbits_of_Solutions is
! 8:
! 9: -- AUXILIAIRIES :
! 10:
! 11: function Equal_upon_Sign
! 12: ( x,y : Standard_Complex_Vectors.Vector; tol : double_float )
! 13: return boolean is
! 14:
! 15: -- DESCRIPTION :
! 16: -- Returns true if for all i : | x(i) +- y(i) | <= tol.
! 17:
! 18: begin
! 19: for i in x'range loop
! 20: if (AbsVal(x(i)-y(i)) > tol) and then (AbsVal(x(i)+y(i)) > tol)
! 21: then return false;
! 22: end if;
! 23: end loop;
! 24: return true;
! 25: end Equal_upon_Sign;
! 26:
! 27: procedure Is_In ( sols : in out Solution_List; s : in Solution;
! 28: sign : in boolean; tol : in double_float;
! 29: isin : out boolean ) is
! 30:
! 31: -- DESCRIPTION :
! 32: -- The output variable isin isin becomes true if one of the solutions
! 33: -- in the list sols equals s. When this is the case, the multiplicity
! 34: -- of the corresponding solution is augmented.
! 35: -- Otherwise, isin equals false on return.
! 36:
! 37: s1 : Solution(s.n);
! 38: temp : Solution_List;
! 39:
! 40: begin
! 41: temp := sols;
! 42: while not Is_Null(temp) loop
! 43: s1 := Head_Of(temp).all;
! 44: if (sign and then Equal_upon_Sign(s.v,s1.v,tol))
! 45: or else (not sign and then Equal(s.v,s1.v,tol))
! 46: then Head_Of(temp).m := Head_Of(temp).m + 1;
! 47: isin := true;
! 48: return;
! 49: end if;
! 50: temp := Tail_Of(temp);
! 51: end loop;
! 52: isin := false;
! 53: end Is_In;
! 54:
! 55: -- CONSTRUCTORS :
! 56:
! 57: procedure Analyze ( l : in List_of_Permutations; sign : in boolean;
! 58: tol : in double_float; sols : in out Solution_List ) is
! 59: begin
! 60: if not Is_Null(sols)
! 61: then declare
! 62: res,res_last,tmpsols : Solution_List;
! 63: n : natural := Head_Of(sols).n;
! 64: s1,s2 : Solution(n);
! 65: lp : Link_to_Permutation;
! 66: tmpl : List_of_Permutations;
! 67: found : boolean;
! 68: begin
! 69: tmpsols := sols;
! 70: while not Is_Null(tmpsols) loop
! 71: s1 := Head_Of(tmpsols).all;
! 72: tmpl := l;
! 73: while not Is_Null(tmpl) loop
! 74: lp := Head_Of(tmpl);
! 75: s2.t := s1.t;
! 76: s2.m := 1;
! 77: s2.v := Permutation(lp.all)*s1.v;
! 78: Is_In(res,s2,sign,tol,found);
! 79: exit when found;
! 80: tmpl := Tail_Of(tmpl);
! 81: end loop;
! 82: if not found
! 83: then Append(res,res_last,s1);
! 84: end if;
! 85: tmpsols := Tail_Of(tmpsols);
! 86: end loop;
! 87: Clear(sols); sols := res;
! 88: end;
! 89: end if;
! 90: end Analyze;
! 91:
! 92: function Permutable ( s : Solution; sign : boolean; tol : double_float;
! 93: sols : Solution_List ) return boolean is
! 94:
! 95: -- DESCRIPTION :
! 96: -- Returns true if the solutions s can be permuted into one of the
! 97: -- solutions of the given list.
! 98:
! 99: tmp : Solution_List := sols;
! 100:
! 101: begin
! 102: while not Is_Null(tmp) loop
! 103: if sign and then Sign_Permutable(s.v,Head_Of(tmp).v,tol)
! 104: then return true;
! 105: elsif Permutable(s.v,Head_Of(tmp).v,tol)
! 106: then return true;
! 107: else tmp := Tail_Of(tmp);
! 108: end if;
! 109: end loop;
! 110: return false;
! 111: end Permutable;
! 112:
! 113: procedure Permutable
! 114: ( s : in Solution; sign : in boolean; tol : in double_float;
! 115: sols : in out Solution_List; isin : out boolean ) is
! 116:
! 117: -- DESCRIPTION :
! 118: -- The output variable isin becomes true when the given solution can be
! 119: -- permuted into one of the solutions of the given list.
! 120: -- In this case, also the multiplicity of the corresponding solutions is
! 121: -- augmented. Otherwise, isin is false on return.
! 122:
! 123: tmp : Solution_List := sols;
! 124: ls1 : Link_to_Solution;
! 125:
! 126: begin
! 127: while not Is_Null(tmp) loop
! 128: ls1 := Head_Of(tmp);
! 129: if (sign and then Sign_Permutable(s.v,ls1.v,tol))
! 130: or else (not sign and then Permutable(s.v,ls1.v,tol))
! 131: then isin := true;
! 132: ls1.m := ls1.m + 1; Set_Head(tmp,ls1);
! 133: return;
! 134: else tmp := Tail_Of(tmp);
! 135: end if;
! 136: end loop;
! 137: isin := false;
! 138: end Permutable;
! 139:
! 140: function Generating ( sols : Solution_List; sign : boolean;
! 141: tol : double_float ) return Solution_List is
! 142:
! 143: tmp,gensols,gensols_last : Solution_List;
! 144:
! 145: begin
! 146: tmp := sols;
! 147: while not Is_Null(tmp) loop
! 148: declare
! 149: ls : Link_to_Solution := Head_Of(tmp);
! 150: isin : boolean;
! 151: begin
! 152: Permutable(ls.all,sign,tol,gensols,isin);
! 153: if not isin
! 154: then Append(gensols,gensols_last,ls.all);
! 155: end if;
! 156: end;
! 157: tmp := Tail_Of(tmp);
! 158: end loop;
! 159: return gensols;
! 160: end Generating;
! 161:
! 162: procedure Orbit_Structure ( s : in Solution; tol : in double_float;
! 163: orbit : in out Permutation;
! 164: nbdiff : out natural ) is
! 165:
! 166: nb : natural := 1;
! 167: found : boolean;
! 168:
! 169: begin
! 170: orbit(1) := nb;
! 171: for i in 2..s.n loop
! 172: found := false;
! 173: for j in 1..(i-1) loop
! 174: if AbsVal(s.v(i) - s.v(j)) < tol
! 175: then orbit(i) := orbit(j);
! 176: found := true;
! 177: exit;
! 178: end if;
! 179: end loop;
! 180: if not found
! 181: then nb := nb + 1;
! 182: orbit(i) := nb;
! 183: end if;
! 184: end loop;
! 185: nbdiff := nb;
! 186: end Orbit_Structure;
! 187:
! 188: function Orbits ( sols : Solution_List; tol : double_float )
! 189: return permutation is
! 190:
! 191: ind : natural;
! 192: tmp : Solution_List := sols;
! 193: n : natural := Head_Of(sols).n;
! 194: o,orb : Permutation(1..n);
! 195: sol : Solution(n);
! 196:
! 197: begin
! 198: orb := (1..n => 0);
! 199: o := orb;
! 200: while not Is_Null(tmp) loop
! 201: sol := Head_Of(tmp).all;
! 202: Orbit_Structure(sol,tol,o,ind);
! 203: orb(ind) := orb(ind) + sol.m;
! 204: tmp := Tail_Of(tmp);
! 205: end loop;
! 206: return orb;
! 207: end Orbits;
! 208:
! 209: procedure Orbits ( grp : in List_of_Permutations; tol : in double_float;
! 210: sols : in out Solution_List;
! 211: lorb : in out List_of_Orbits ) is
! 212:
! 213: tmp : Solution_List;
! 214: n : natural := Head_Of(sols).n;
! 215:
! 216: begin
! 217: Analyze(grp,false,tol,sols);
! 218: tmp := sols;
! 219: while not Is_Null(tmp) loop
! 220: declare
! 221: sol : Solution(n) := Head_Of(tmp).all;
! 222: orbi : Orbit(n);
! 223: begin
! 224: Orbit_Structure(sol,tol,orbi.orb,orbi.nbdiff);
! 225: orbi.nbgen := 1;
! 226: orbi.nbsols := sol.m;
! 227: declare
! 228: tmp2 : List_of_Orbits := lorb;
! 229: found : boolean := false;
! 230: begin
! 231: while not Is_Null(tmp2) loop
! 232: declare
! 233: orbi2 : Orbit(n) := Head_Of(tmp2).all;
! 234: lorbi2 : Link_to_Orbit := Head_Of(tmp2);
! 235: begin
! 236: if orbi2.nbdiff = orbi.nbdiff
! 237: and then Same_Orbit(orbi.orb,orbi2.orb)
! 238: then lorbi2.nbgen := orbi2.nbgen + 1;
! 239: lorbi2.nbsols := orbi2.nbsols + sol.m;
! 240: found := true;
! 241: exit;
! 242: else tmp2 := Tail_Of(tmp2);
! 243: end if;
! 244: end;
! 245: end loop;
! 246: if not found
! 247: then declare
! 248: liorbi : Link_to_Orbit := new Orbit'(orbi);
! 249: begin
! 250: Construct(liorbi,lorb);
! 251: end;
! 252: end if;
! 253: end;
! 254: end;
! 255: tmp := Tail_Of(tmp);
! 256: end loop;
! 257: end Orbits;
! 258:
! 259: -- SELECTOR :
! 260:
! 261: function Same_Orbit ( orb1,orb2 : Permutation ) return boolean is
! 262:
! 263: f1,f2 : Permutation(1..orb1'last);
! 264: found : boolean;
! 265:
! 266: begin
! 267: -- construct frequencies :
! 268: f1 := (f1'range => 0);
! 269: for i in orb1'range loop
! 270: f1(orb1(i)) := f1(orb1(i)) + 1;
! 271: end loop;
! 272: f2 := (f2'range => 0);
! 273: for i in orb2'range loop
! 274: f2(orb2(i)) := f2(orb2(i)) + 1;
! 275: end loop;
! 276: -- compare the frequencies :
! 277: for i in f1'range loop
! 278: found := false;
! 279: for j in f2'range loop
! 280: if f1(i) = f2(j)
! 281: then found := true;
! 282: exit;
! 283: end if;
! 284: end loop;
! 285: if not found
! 286: then return false;
! 287: end if;
! 288: end loop;
! 289: for i in f2'range loop
! 290: found := false;
! 291: for j in f1'range loop
! 292: if f2(i) = f1(j)
! 293: then found := true;
! 294: exit;
! 295: end if;
! 296: end loop;
! 297: if not found
! 298: then return false;
! 299: end if;
! 300: end loop;
! 301: return true;
! 302: end Same_Orbit;
! 303:
! 304: -- DESTRUCTOR :
! 305:
! 306: procedure Clear ( lorb : in out List_of_Orbits ) is
! 307:
! 308: procedure free is new unchecked_deallocation(Orbit,Link_to_Orbit);
! 309:
! 310: tmp : List_of_Orbits := lorb;
! 311: begin
! 312: while not Is_Null(tmp) loop
! 313: declare
! 314: lior : Link_to_Orbit := Head_Of(tmp);
! 315: begin
! 316: free(lior);
! 317: end;
! 318: tmp := Tail_Of(tmp);
! 319: end loop;
! 320: Lists_of_Orbits.Clear(Lists_of_Orbits.List(lorb));
! 321: end Clear;
! 322:
! 323: end Orbits_of_Solutions;
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>