[BACK]Return to orbits_of_solutions.adb CVS log [TXT][DIR] Up to [local] / OpenXM_contrib / PHC / Ada / Root_Counts / Symmetry

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>