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

Annotation of OpenXM_contrib/PHC/Ada/Root_Counts/Stalift/mixed_volume_computation.adb, Revision 1.1

1.1     ! maekawa     1: with integer_io;                         use integer_io;
        !             2: with Standard_Integer_Vectors_io;        use Standard_Integer_Vectors_io;
        !             3: with Arrays_of_Integer_Vector_Lists_io;  use Arrays_of_Integer_Vector_Lists_io;
        !             4:
        !             5: with Standard_Floating_Vectors;
        !             6: with Standard_Integer_Matrices;          use Standard_Integer_Matrices;
        !             7: with Standard_Integer_Linear_Solvers;    use Standard_Integer_Linear_Solvers;
        !             8: with Lists_of_Integer_Vectors;           use Lists_of_Integer_Vectors;
        !             9: with Transforming_Integer_Vector_Lists;  use Transforming_Integer_Vector_Lists;
        !            10: with Integer_Lifting_Utilities;          use Integer_Lifting_Utilities;
        !            11: with Integer_Mixed_Subdivisions_io;      use Integer_Mixed_Subdivisions_io;
        !            12: with Mixed_Coherent_Subdivisions;        use Mixed_Coherent_Subdivisions;
        !            13:
        !            14: package body Mixed_Volume_Computation is
        !            15:
        !            16: -- AUXILIAIRY OUTPUT ROUTINES :
        !            17:
        !            18:   procedure put ( file : in file_type; points : in Array_of_Lists;
        !            19:                   n : in natural; mix : in Vector;
        !            20:                   mixsub : in out Mixed_Subdivision; mv : out natural ) is
        !            21:   begin
        !            22:     new_line(file);
        !            23:     put_line(file,"THE LIFTED SUPPORTS :");
        !            24:     new_line(file);
        !            25:     put(file,points);
        !            26:     new_line(file);
        !            27:     put_line(file,"THE MIXED SUBDIVISION :");
        !            28:     new_line(file);
        !            29:     put(file,n,mix,mixsub,mv);
        !            30:   end put;
        !            31:
        !            32:   procedure Sort ( supports : in out Array_of_Lists; k,nb,n : in natural;
        !            33:                    mxt,perm : in out Vector ) is
        !            34:
        !            35:   -- DESCRIPTION :
        !            36:   --   Auxiliary operation for Compute_Mixture.
        !            37:   --   Compares the kth support with the following supports.
        !            38:   --   Already nb different supports have been found.
        !            39:
        !            40:   begin
        !            41:     for l in (k+1)..n loop
        !            42:       if Equal(Supports(k),Supports(l))
        !            43:        then if l /= k + mxt(nb)
        !            44:              then declare
        !            45:                     pos : natural := k + mxt(nb);
        !            46:                     tmpdl : List := supports(l);
        !            47:                     tmppos : natural;
        !            48:                   begin
        !            49:                     supports(l) := supports(pos);
        !            50:                     supports(pos) := tmpdl;
        !            51:                     tmppos := perm(l);
        !            52:                     perm(l) := perm(pos);
        !            53:                     perm(pos) := tmppos;
        !            54:                   end;
        !            55:             end if;
        !            56:             mxt(nb) := mxt(nb) + 1;
        !            57:       end if;
        !            58:     end loop;
        !            59:   end Sort;
        !            60:
        !            61: -- TARGET ROUTINES :
        !            62:
        !            63:   procedure Compute_Mixture ( supports : in out Array_of_Lists;
        !            64:                               mix,perms : out Link_to_Vector ) is
        !            65:
        !            66:     n : constant natural := supports'last;
        !            67:     cnt : natural := 0;              -- counts the number of different supports
        !            68:     mxt : Vector(supports'range)     -- counts the number of occurrencies
        !            69:         := (supports'range => 1);
        !            70:     perm : Link_to_Vector            -- keeps track of the permutations
        !            71:          := new Standard_Integer_Vectors.Vector(supports'range);
        !            72:     index : natural := supports'first;
        !            73:
        !            74:   begin
        !            75:     for k in perm'range loop
        !            76:       perm(k) := k;
        !            77:     end loop;
        !            78:     while index <= supports'last loop
        !            79:       cnt := cnt + 1;
        !            80:       Sort(supports,index,cnt,n,mxt,perm.all);
        !            81:       index := index + mxt(cnt);
        !            82:     end loop;
        !            83:     mix := new Standard_Integer_Vectors.Vector'(mxt(mxt'first..cnt));
        !            84:     perms := perm;
        !            85:   end Compute_Mixture;
        !            86:
        !            87:   function Compute_Index ( k : natural; mix : Vector ) return natural is
        !            88:
        !            89:   -- DESCRIPTION :
        !            90:   --   Returns the index of k w.r.t. to the type of mixture.
        !            91:
        !            92:     index : natural := mix(mix'first);
        !            93:
        !            94:   begin
        !            95:     if k <= index
        !            96:      then return mix'first;
        !            97:      else for l in (mix'first+1)..mix'last loop
        !            98:             index := index + mix(l);
        !            99:             if k <= index
        !           100:              then return l;
        !           101:             end if;
        !           102:           end loop;
        !           103:           return mix'last;
        !           104:     end if;
        !           105:   end Compute_Index;
        !           106:
        !           107:   function Compute_Permutation ( n : natural; mix : Vector;
        !           108:                                  supports : Array_of_Lists )
        !           109:                                return Link_to_Vector is
        !           110:
        !           111:     perms : Link_to_Vector := new Vector(1..n);
        !           112:
        !           113:   begin
        !           114:     for k in perms'range loop
        !           115:       perms(k) := k;
        !           116:     end loop;
        !           117:     return perms;
        !           118:   end Compute_Permutation;
        !           119:
        !           120:   function Permute ( p : Poly_Sys; perm : Link_to_Vector ) return Poly_Sys is
        !           121:
        !           122:     res : Poly_Sys(p'range);
        !           123:
        !           124:   begin
        !           125:     for k in p'range loop
        !           126:       res(k) := p(perm(k));
        !           127:     end loop;
        !           128:     return res;
        !           129:   end Permute;
        !           130:
        !           131:   function Permute ( supports : Array_of_Lists; perm : Link_to_Vector )
        !           132:                    return Array_of_Lists is
        !           133:
        !           134:     res : Array_of_Lists(supports'range);
        !           135:
        !           136:   begin
        !           137:     for k in supports'range loop
        !           138:       res(k) := supports(perm(k));
        !           139:     end loop;
        !           140:     return res;
        !           141:   end Permute;
        !           142:
        !           143:   function Typed_Lists ( mix : Vector; points : Array_of_Lists )
        !           144:                        return Array_of_Lists is
        !           145:
        !           146:     res : Array_of_Lists(mix'range);
        !           147:     ind : natural := res'first;
        !           148:
        !           149:   begin
        !           150:     for i in mix'range loop
        !           151:       res(i) := points(ind);
        !           152:       ind := ind + mix(i);
        !           153:     end loop;
        !           154:     return res;
        !           155:   end Typed_Lists;
        !           156:
        !           157: -- MIXED VOLUME COMPUTATIONS BASED ON SUBDIVISIONS :
        !           158:
        !           159: -- AUXILIARIES :
        !           160:
        !           161:   function Is_Fine ( mix : Vector; mic : Mixed_Cell ) return boolean is
        !           162:
        !           163:   -- DESCRIPTION :
        !           164:   --   Returns true if the mixed volume can be computed by a determinant.
        !           165:
        !           166:     fine : boolean := true;
        !           167:
        !           168:   begin
        !           169:     for k in mic.pts'range loop
        !           170:       fine := (Length_Of(mic.pts(k)) = mix(k) + 1);
        !           171:       exit when not fine;
        !           172:     end loop;
        !           173:     return fine;
        !           174:   end Is_Fine;
        !           175:
        !           176:   function Reduced_Supports ( n : natural; mix : Vector; mic : Mixed_Cell )
        !           177:                             return Array_of_Lists is
        !           178:
        !           179:   -- DESCRIPTION :
        !           180:   --   Returns the supports of the cell without the lifting values.
        !           181:
        !           182:     res : Array_of_Lists(1..n);
        !           183:     cnt : natural := 1;
        !           184:
        !           185:   begin
        !           186:     for k in mic.pts'range loop
        !           187:       res(cnt) := Reduce(mic.pts(k),n+1);
        !           188:       for l in 1..mix(k)-1 loop
        !           189:         Copy(res(cnt),res(cnt+l));
        !           190:       end loop;
        !           191:       cnt := cnt + mix(k);
        !           192:     end loop;
        !           193:     return res;
        !           194:   end Reduced_Supports;
        !           195:
        !           196:   function Fine_Mixed_Volume ( n : natural; mix : Vector; mic : Mixed_Cell )
        !           197:                              return natural is
        !           198:
        !           199:   -- DESCRIPTION :
        !           200:   --   Computes the mixed volume for a cell that is fine mixed.
        !           201:
        !           202:   -- REQUIRED : Fine(mix,mic).
        !           203:
        !           204:     res,count : natural;
        !           205:     mat : matrix(1..n,1..n);
        !           206:     detmat : integer;
        !           207:     tmp : List;
        !           208:     sh,pt : Link_to_Vector;
        !           209:
        !           210:   begin
        !           211:     count := 1;
        !           212:     for k in mic.pts'range loop
        !           213:       sh := Head_Of(mic.pts(k));
        !           214:       tmp := Tail_Of(mic.pts(k));
        !           215:       while not Is_Null(tmp) loop
        !           216:         pt := Head_Of(tmp);
        !           217:         for j in 1..n loop
        !           218:           mat(count,j) := pt(j) - sh(j);
        !           219:         end loop;
        !           220:         tmp := Tail_Of(tmp);
        !           221:         count := count + 1;
        !           222:       end loop;
        !           223:     end loop;
        !           224:     detmat := Det(mat);
        !           225:     if detmat >= 0
        !           226:      then res := detmat;
        !           227:      else res := -detmat;
        !           228:     end if;
        !           229:     return res;
        !           230:   end Fine_Mixed_Volume;
        !           231:
        !           232:   function Mixed_Volume ( n : natural; mix : Vector;
        !           233:                          mic : Mixed_Cell ) return natural is
        !           234:
        !           235:   -- ALGORITHM :
        !           236:   --   First check if the cell has a refinement, if so, then use it,
        !           237:   --   if not, then check if the cell is fine mixed.
        !           238:   --   If the cell is fine mixed, only a determinant needs to be computed,
        !           239:   --   otherwise the cell will be refined.
        !           240:
        !           241:     res : natural;
        !           242:
        !           243:   begin
        !           244:     if (mic.sub /= null) and then not Is_Null(mic.sub.all)
        !           245:      then res := Mixed_Volume_Computation.Mixed_Volume(n,mix,mic.sub.all);
        !           246:      elsif Is_Fine(mix,mic)
        !           247:          then res := Fine_Mixed_Volume(n,mix,mic);
        !           248:          else declare
        !           249:                 rcell : Array_of_Lists(1..n) := Reduced_Supports(n,mix,mic);
        !           250:               begin
        !           251:                 res := Mixed_Volume_Computation.Mixed_Volume(n,rcell);
        !           252:                 Deep_Clear(rcell);
        !           253:               end;
        !           254:     end if;
        !           255:     return res;
        !           256:   end Mixed_Volume;
        !           257:
        !           258:   function Mixed_Volume ( n : natural; mix : Vector;
        !           259:                          mixsub : Mixed_Subdivision ) return natural is
        !           260:
        !           261:     res : natural := 0;
        !           262:     tmp : Mixed_Subdivision := mixsub;
        !           263:
        !           264:   begin
        !           265:     while not Is_Null(tmp) loop
        !           266:       res := res + Mixed_Volume(n,mix,Head_Of(tmp));
        !           267:       tmp := Tail_Of(tmp);
        !           268:     end loop;
        !           269:     return res;
        !           270:   end Mixed_Volume;
        !           271:
        !           272:   procedure Mixed_Volume ( n : in natural; mix : in Vector;
        !           273:                            mic : in out Mixed_Cell; mv : out natural ) is
        !           274:   begin
        !           275:     if (mic.sub /= null) and then not Is_Null(mic.sub.all)
        !           276:      then mv := Mixed_Volume_Computation.Mixed_Volume(n,mix,mic.sub.all);
        !           277:      elsif Is_Fine(mix,mic)
        !           278:          then mv := Fine_Mixed_Volume(n,mix,mic);
        !           279:          else -- NOTE : keep the same type of mixture!
        !           280:            declare
        !           281:              rcell : Array_of_Lists(1..n) := Reduced_Supports(n,mix,mic);
        !           282:              lifted : Array_of_Lists(mix'range);
        !           283:              mixsub : Mixed_Subdivision;
        !           284:            begin
        !           285:              Mixed_Volume_Computation.Mixed_Volume
        !           286:                (n,mix,rcell,lifted,mixsub,mv);
        !           287:              mic.sub := new Mixed_Subdivision'(mixsub);
        !           288:              Deep_Clear(rcell);  Deep_Clear(lifted);
        !           289:            end;
        !           290:     end if;
        !           291:   end Mixed_Volume;
        !           292:
        !           293:   procedure Mixed_Volume ( n : in natural; mix : in Vector;
        !           294:                            mixsub : in out Mixed_Subdivision;
        !           295:                            mv : out natural ) is
        !           296:
        !           297:     res : natural := 0;
        !           298:     tmp : Mixed_Subdivision := mixsub;
        !           299:
        !           300:   begin
        !           301:     while not Is_Null(tmp) loop
        !           302:       declare
        !           303:         mic : Mixed_Cell := Head_Of(tmp);
        !           304:         mmv : natural;
        !           305:       begin
        !           306:         Mixed_Volume(n,mix,mic,mmv);
        !           307:         Set_Head(tmp,mic);
        !           308:         res := res + mmv;
        !           309:       end;
        !           310:       tmp := Tail_Of(tmp);
        !           311:     end loop;
        !           312:     mv := res;
        !           313:   end Mixed_Volume;
        !           314:
        !           315: -- MIXED VOLUME COMPUTATIONS BASED ON SUPPORTS :
        !           316:
        !           317:   function Mixed_Volume ( n : natural; supports : Array_of_Lists )
        !           318:                        return natural is
        !           319:     mv : natural;
        !           320:     mix,perm : Link_to_Vector;
        !           321:     permsupp : Array_of_Lists(supports'range);
        !           322:
        !           323:   begin
        !           324:     Copy(supports,permsupp);
        !           325:     Compute_Mixture(permsupp,mix,perm);
        !           326:     mv := Mixed_Volume(n,mix.all,permsupp);
        !           327:     Clear(mix); Clear(perm); Deep_Clear(permsupp);
        !           328:     return mv;
        !           329:   end Mixed_Volume;
        !           330:
        !           331:   function Mixed_Volume ( file : file_type; n : natural;
        !           332:                           supports : Array_of_Lists ) return natural is
        !           333:     mv : natural;
        !           334:     mix,perm : Link_to_Vector;
        !           335:     permsupp : Array_of_Lists(supports'range);
        !           336:
        !           337:   begin
        !           338:     Copy(supports,permsupp);
        !           339:     Compute_Mixture(permsupp,mix,perm);
        !           340:     mv := Mixed_Volume(file,n,mix.all,permsupp);
        !           341:     Clear(mix); Clear(perm); Deep_Clear(permsupp);
        !           342:     return mv;
        !           343:   end Mixed_Volume;
        !           344:
        !           345:   function Mixed_Volume ( n : natural; mix : Vector;
        !           346:                           supports : Array_of_Lists ) return natural is
        !           347:
        !           348:     res : natural;
        !           349:     mixsub : Mixed_Subdivision;
        !           350:     lifted : Array_of_Lists(mix'range);
        !           351:
        !           352:   begin
        !           353:     Mixed_Volume_Computation.Mixed_Volume(n,mix,supports,lifted,mixsub,res);
        !           354:     Deep_Clear(lifted);
        !           355:     Shallow_Clear(mixsub);
        !           356:     return res;
        !           357:   end Mixed_Volume;
        !           358:
        !           359:   procedure Mixed_Volume
        !           360:               ( n : in natural; mix : in Vector; supports : in Array_of_Lists;
        !           361:                 lifted : out Array_of_Lists;
        !           362:                 mixsub : out Mixed_Subdivision; mv : out natural ) is
        !           363:
        !           364:     low : constant Vector := (mix'range => 0);
        !           365:     upp : constant Vector := Adaptive_Lifting(supports);
        !           366:     nbsucc,nbfail : Standard_Floating_Vectors.Vector(mix'range)
        !           367:                   := (mix'range => 0.0);
        !           368:     liftsupp : Array_of_Lists(mix'range);
        !           369:     sub : Mixed_Subdivision;
        !           370:
        !           371:   begin
        !           372:     Mixed_Coherent_Subdivision(n,mix,supports,false,low,upp,liftsupp,
        !           373:                                nbsucc,nbfail,sub);
        !           374:     Mixed_Volume(n,mix,sub,mv);
        !           375:     lifted := liftsupp;
        !           376:     mixsub := sub;
        !           377:   end Mixed_Volume;
        !           378:
        !           379:   function Mixed_Volume ( file : file_type; n : natural; mix : Vector;
        !           380:                           supports : Array_of_Lists ) return natural is
        !           381:
        !           382:     res : natural;
        !           383:     mixsub : Mixed_Subdivision;
        !           384:     lifted : Array_of_Lists(mix'range);
        !           385:
        !           386:   begin
        !           387:     Mixed_Volume_Computation.Mixed_Volume
        !           388:       (file,n,mix,supports,lifted,mixsub,res);
        !           389:     Deep_Clear(lifted);
        !           390:     Shallow_Clear(mixsub);
        !           391:     return res;
        !           392:   end Mixed_Volume;
        !           393:
        !           394:   procedure Mixed_Volume
        !           395:               ( file : in file_type; n : in natural;
        !           396:                 mix : in Vector; supports : in Array_of_Lists;
        !           397:                 lifted : out Array_of_Lists;
        !           398:                 mixsub : out Mixed_Subdivision; mv : out natural ) is
        !           399:
        !           400:     low : constant Vector := (mix'range => 0);
        !           401:     upp : constant Vector := Adaptive_Lifting(supports);
        !           402:     sub : Mixed_Subdivision;
        !           403:     nbsucc,nbfail : Standard_Floating_Vectors.Vector(mix'range)
        !           404:                   := (mix'range => 0.0);
        !           405:     liftsupp : Array_of_Lists(mix'range);
        !           406:
        !           407:   begin
        !           408:     Mixed_Coherent_Subdivision
        !           409:       (n,mix,supports,false,low,upp,liftsupp,nbsucc,nbfail,sub);
        !           410:     lifted := liftsupp;
        !           411:     mixsub := sub;
        !           412:     put(file,liftsupp,n,mix,sub,mv);
        !           413:   end Mixed_Volume;
        !           414:
        !           415: end Mixed_Volume_Computation;

FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>