progp

Best algorithm for LZW ..C language

      if not (edge in [hin,hout]) then
      begin
         OK:=false;
	 writeln('**hi edge value incorrect');
      end;
   end;
   
   CheckHi:=OK;
end;{CheckHi}
  
function CheckLo(X:Sreal):boolean;
var OK:boolean;
begin
   OK:=true;
   with X do
   begin
      case cardinality of
      infinite:
         if (exp=Maxexp)and(mantissa=Mininf) then
	 else writeln('**Invalid lo infinity');
      finite:
      begin
         if (mantissa=Maxinf) or (mantissa=Mininf) then
	 begin OK:=false; writeln('**Invalid finite value - hi');
	 end;
	 
         if mantissa = 0 then
	    if (exp=0) then 
	    else 
	    begin OK:=false; writeln('**Invalid zero - lo')
	    end
	 else
	 begin
	       if (mantissa > 0) then
	          if mantissa >= (Maxinf div 10) then{OK}
		  else 
		  begin OK:=false; writeln('**Incorrect normalization - lo') 
		  end
	       else{mantissa<0}
	          if mantissa > (Mininf div 10) then
		  begin OK:=false; writeln('**Incorrect normalization - lo') 
		  end;
	 end;
      end;
      end;{case}

      if not (edge in [lin,lout]) then
      begin
         OK:=false;
	 writeln('**lo edge value incorrect');
      end;
   end;
   
   CheckLo:=OK;
end;{CheckLo}
  
function CheckInt(I:Int):boolean;
var OK:boolean;
begin
   OK:=CheckHi(I.hi) and CheckLo(I.lo);
   if gtS(I.lo,I.hi) then
   begin
      OK:=false;
      writeln('**Limits out of order');
   end;

   if not OK then 
   begin writeln('**Error in Check'); DumpInt(I);
   end;
   
   CheckInt:=OK;
end;

procedure DumpMem(var DCurr:DataMem);
var tD:Ptr; tL:Loc;
begin
   with DCurr do
   begin
        writeln('LastHalve:',LastHalve:0);
	
    	for tL:= 1 to End do
	   writeln(tL:3,S[tL]:2,RHalve[tL]);
	writeln;
	
	for tD:= 1 to MaxDMem do 
	begin 
	   write(tD:5);
	   DumpInt(D[tD]);
	   assert(CheckInt(D[tD]));
	   writeln;
	end;
	writeln;
   end;
end;{DumpMem}

procedure WriteMem(var DCurr:DataMem);
var tD:Ptr; 
begin
   with DCurr do
   begin
	for tD:= 1 to MaxDMem do 
	if (DF.PF[tD] > PNull) or (Debug > activity) then
	begin 
	   write(tD:5);
	   WriteInt(D[tD]);
	   writeln;
	end;
	writeln;
   end;
end;{WriteMem}

procedure OuterExec
(PC:Loc0;DCurr:DataMem;Change:boolean;First:State;
 var OldCounter:Positive;Level:Positive);

var Counter:Positive;
    Fail,AllPoints,LocalChange:boolean;

procedure NewOuter(F:State);
begin OuterExec(PC,DCurr,Change,F,Counter,Level+1);
end;

{!!}
procedure execprint(PC:Loc; L:Ptr; R0:Int);
begin
   DF.PF[L]:=PSoln;
   writeln;
   write(PC:3,L:5);
   WriteInt(R0);
   writeln;
end;

procedure execpr(var Sr:State; L:Ptr);
begin
   Sr:=-1; DF.PF[L]:=PPrint;
end;{execpr}

procedure exectr(var Sr:State; L:Ptr);
begin
   Sr:=-1; DF.PF[L]:=PTrace;
end;{exectr}

procedure execsoln(var Sr:State; L:Ptr);
begin
   Sr:=-1; DF.PF[L]:=PSoln;
end;{execsoln}

procedure execreadr(var Sr:State;var R0:Int);
begin
   writeln;
   write('<<');
   ReadInt(R0);
   Sr:=-1;
end;

function GetReal(E,M:integer):real;
{convert E-exponent,M-mantissa into genuine Pascal real number}
var x:real;
begin
   x:=M/Maxinf;
   while E>0 do begin x:=x*10; E:=E-1; end; 
   while E<0 do begin x:=x/10; E:=E+1; end;
   GetReal:=x; 
end;{GetReal} 
   
procedure Ratio(Lo,Hi:Sreal;var ERat,MRat:integer); 
{compute ratio of Hi to Lo in exponent mantissa form}
begin
   if Lo.mantissa=0 then
   begin{treat zero as if smallest possible positive number}
      ERat:=Hi.exp-Minexp;
      MRat:=Hi.mantissa*10;
   end else
   if Hi.mantissa=0 then
   begin{treat zero as if smallest possible negative number}
      ERat:=Minexp-Lo.exp;
      MRat:=Lo.mantissa*10;
   end
   else begin
      ERat:=Hi.exp-Lo.exp;
      MRat:=(Hi.mantissa*Maxinf) div Lo.mantissa;
   end;
end;{Ratio}
    
function Adjacent(X:Int):boolean;
{are hi and lo bounds adjacent points}
begin
   with X do
   if (hi.mantissa=0) or (lo.mantissa=0) then
      Adjacent:=
         ((hi.mantissa=0)and(lo.mantissa=Mininf div 10)and(lo.exp=Minexp)) or
         ((lo.mantissa=0)and(hi.mantissa=Maxinf div 10)and(hi.exp=Minexp)) 
   else
      Adjacent:=
         ((lo.exp=hi.exp)and(lo.mantissa+1=hi.mantissa)) or
	 ((hi.exp=lo.exp+1)and(hi.mantissa=(lo.mantissa div 10)+1)) or
	 ((hi.exp=lo.exp-1)and((hi.mantissa div 10)-1=lo.mantissa));
end;{Adjacent}

procedure exechalve
   (var PC:Loc0;var Sr:State;var R0:Int;var OK:boolean;var Change:boolean);
{Reduce range of R0 (suceeds twice for two 'halves')}

var EDiff,MDiff,ERat,MRat,MidE,MidM,M0,M1,HiM,HiE:integer;
    Dummy:boolean;
    Mid:Sreal;
    R,D:real;
    OldPC:Loc;
   
    procedure AtEnd;{What to do afer a successful halve}
    begin  
	DCurr.LastHalve:=PC; PC:=0; Sr:=0;
    end;

begin{exechalve}
OldPC:=PC;
with R0 do
 begin
      if DCurr.LastHalve >= PC then {not our turn yet} else
      if (lo.mantissa = hi.mantissa) and (lo.exp=hi.exp) and
         (lo.edge=lin) and (hi.edge=hin)
      then {single point cant be divided} Sr:=-1 
      else
      if Adjacent(R0) and 
	 (((lo.edge=lout) and (hi.edge=hout)) or
	  ((lo.cardinality=infinite)and(hi.edge=hout)) or 
	  ((hi.cardinality=infinite)and(lo.edge=lout))
	 )
      then Sr:=-1
      else
      if Sr=0 then
         begin
            AlignUp(hi.exp,hi.mantissa,lo.exp,-lo.mantissa,EDiff,M0,M1,Dummy);
            MDiff:=M0+M1;
	    D:=GetReal(EDiff,MDiff);
	    if D < DCurr.RHalve[PC]*HalveThreshold 
	    then {already narrowed enough dont bother} 
	       Change:=true {otherwise can terminate too early}
	    else begin
	       DCurr.RHalve[PC]:=D;
               NewOuter(1);
               NewOuter(2);
	       OK:=false;{fail after both alternatives tried}
	    end
         end
   else{Sr=1,2}
   if Adjacent(R0)
   then begin{two adjacent points - needs special care}
      if (hi.edge=hin) and (hi.cardinality=finite) then
      begin
         case Sr of
	 1:hi.edge:=hout;
	 2:begin lo:=hi; lo.edge:=lin;
	   end;
	 end;

         AtEnd;
      end
      else if (lo.edge=lin) and (lo.cardinality=finite) then 
      begin
         case Sr of
	 1:begin hi:=lo; hi.edge:=hin;
	   end;
	 2:lo.edge:=lout;
	 end;
         AtEnd;
      end else {cant be narrowed} Sr:=-1;

   end{adjacent} else
   begin
      if (lo.mantissa < 0) and (hi.mantissa > 0) then
      begin
         MidM:=0; MidE:=0;
      end else
      begin
         Ratio(lo,hi,ERat,MRat);
         if MRat < 0 then MRat:=-MRat;
         R:=GetReal(ERat,MRat);
(*writeln(ERat,MRat,R);*)
         AlignUp(hi.exp,hi.mantissa,lo.exp,-lo.mantissa,EDiff,M0,M1,Dummy);
         MDiff:=M0+M1;
         if (R > 4) or (R < 0.25) then
         begin{divide hi by sqrt of ratio to get midpoint}
	    if hi.mantissa = 0 
	    then begin HiM:=Mininf div 10; HiE:=Minexp;
	    end else
	    begin HiM:=hi.mantissa; HiE:=hi.exp;
	    end;
	    if ERat < 0 then MidE:=HiE-((ERat-1) div 2)
	 	        else MidE:=HiE-(ERat div 2);
	    if odd(ERat) 
	       then MidM:=trunc(HiM*(Maxinf div 100)/sqrt(MRat*10))
	       else MidM:=trunc(HiM*(Maxinf div 100)/sqrt(MRat));
(*writeln(MidE,MidM);*)
         end else
         begin{take (hi+lo)/2 as midpoint}
            MidM:=MDiff div 2 - M1;
	    MidE:=EDiff;
         end;
      end;
      if MidM >= 0 then NormalizeDn(MidE,MidM,Mid,Dummy)
      		   else NormalizeUp(MidE,MidM,Mid,Dummy);
      case Sr of
      1:begin hi:=Mid; hi.edge:=hout;
        end;
      2:begin lo:=Mid; lo.edge:=lin;
        end;
      end;
(*DumpInt(R0);writeln;*)
      AtEnd;
   end;{if Sr}

   AlignUp(hi.exp,hi.mantissa,lo.exp,-lo.mantissa,EDiff,M0,M1,Dummy);
   MDiff:=M0+M1;
   DCurr.RHalve[OldPC]:=GetReal(EDiff,MDiff);

 end;{with}
end;{exechalve}

procedure exechalves
   (var PC:Loc0;var Sr:State;var R0:Int;var OK:boolean;var Change:boolean);
{Reduce range of R0 (suceeds twice for two 'halves')}
{Simple version thats averages exponents}

var EDiff,MDiff,ERat,MRat,MidE,MidM,M0,M1,HiM,HiE:integer;
    Dummy:boolean;
    Mid:Sreal;
    R,D:real;
    OldPC:Loc;
   
    procedure AtEnd;{What to do afer a successful halve}
    begin  
	DCurr.LastHalve:=PC; PC:=0; Sr:=0;
    end;

    procedure Average(Lo,Hi:Sreal;var Exp:integer);
    {compute average of exponents allowing for zero}
    {infinities happen to work because of representation}
    var Le,He:integer;
    begin
       if Lo.mantissa = 0 then
	  Le := Minexp
       else
	  Le := Lo.exp;
       if Hi.mantissa = 0 then
	  He := Minexp
       else
	  He := Hi.exp;
       Exp:= (He + Le - 2*Minexp) div 2 + Minexp;
writeln(Exp,Hi.exp,Lo.exp,Minexp,He,Le);
    end;{Average}

begin{exechalves}
OldPC:=PC;
with R0 do
 begin
      if DCurr.LastHalve >= PC then {not our turn yet} else
      if (lo.mantissa = hi.mantissa) and (lo.exp=hi.exp) and
         (lo.edge=lin) and (hi.edge=hin)
      then {single point cant be divided} Sr:=-1 
      else
      if Adjacent(R0) and 
	 (((lo.edge=lout) and (hi.edge=hout)) or
	  ((lo.cardinality=infinite)and(hi.edge=hout)) or 
	  ((hi.cardinality=infinite)and(lo.edge=lout))
	 )
      then Sr:=-1
      else
      if Sr=0 then
         begin
            AlignUp(hi.exp,hi.mantissa,lo.exp,-lo.mantissa,EDiff,M0,M1,Dummy);
            MDiff:=M0+M1;
	    D:=GetReal(EDiff,MDiff);
	    if D < DCurr.RHalve[PC]*HalveThreshold 
	    then {already narrowed enough dont bother} 
	       Change:=true {otherwise can terminate too early}
	    else begin
	       DCurr.RHalve[PC]:=D;
               NewOuter(1);
               NewOuter(2);
	       OK:=false;{fail after both alternatives tried}
	    end
         end
   else{Sr=1,2}
   if Adjacent(R0)
   then begin{two adjacent points - needs special care}
      if (hi.edge=hin) and (hi.cardinality=finite) then
      begin
         case Sr of
	 1:hi.edge:=hout;
	 2:begin lo:=hi; lo.edge:=lin;
	   end;
	 end;

         AtEnd;
      end
      else if (lo.edge=lin) and (lo.cardinality=finite) then 
      begin
         case Sr of
	 1:begin hi:=lo; hi.edge:=hin;
	   end;
	 2:lo.edge:=lout;
	 end;
         AtEnd;
      end else {cant be narrowed} Sr:=-1;

   end{adjacent} else
   begin
      if (lo.mantissa < 0) and (hi.mantissa > 0) then
      begin
         MidM:=0; MidE:=0;
      end else
      begin
         Ratio(lo,hi,ERat,MRat);
(*writeln(ERat,MRat,R);*)
         AlignUp(hi.exp,hi.mantissa,lo.exp,-lo.mantissa,EDiff,M0,M1,Dummy);
         MDiff:=M0+M1;
         if (ERat > 1) or (ERat < -1) then
         begin{Average exponents}
	    if hi.mantissa <= 0 
	    then begin MidM:= -Splitman;
	    end else
	    begin MidM:= Splitman; assert(lo.mantissa >= 0);
	    end;
            Average(lo,hi,MidE);
         end else
         begin{take (hi+lo)/2 as midpoint}
            MidM:=MDiff div 2 - M1;
	    MidE:=EDiff;
         end;
      end;
      if MidM >= 0 then NormalizeDn(MidE,MidM,Mid,Dummy)
      		   else NormalizeUp(MidE,MidM,Mid,Dummy);
      case Sr of
      1:begin lo:=Mid; lo.edge:=lin;
        end;
      2:begin hi:=Mid; hi.edge:=hout;
        end;
      end;
(*DumpInt(R0);writeln;*)
      AtEnd;
   end;{if Sr}

   AlignUp(hi.exp,hi.mantissa,lo.exp,-lo.mantissa,EDiff,M0,M1,Dummy);
   MDiff:=M0+M1;
   DCurr.RHalve[OldPC]:=GetReal(EDiff,MDiff);

 end;{with}
end;{exechalves}

procedure execlinh
   (var PC:Loc0;var Sr:State;var R0:Int;var OK:boolean;var Change:boolean);
{Reduce range of R0 (suceeds twice for two 'halves')}

var EDiff,MDiff,MidE,MidM,M0,M1:integer;
    Dummy:boolean;
    Mid:Sreal;
    D:real;
    OldPC:Loc;
   
    procedure AtEnd;{What to do afer a successful halve}
    begin  
	DCurr.LastHalve:=PC; PC:=0; Sr:=0;
    end;

begin{execlinh}
OldPC:=PC;
with R0 do
 begin
      if DCurr.LastHalve >= PC then {not our turn yet} else
      if (lo.mantissa = hi.mantissa) and (lo.exp=hi.exp) and
         (lo.edge=lin) and (hi.edge=hin)
      then {single point cant be divided} Sr:=-1 
      else
      if Adjacent(R0) and 
	 (((lo.edge=lout) and (hi.edge=hout)) or
	  ((lo.cardinality=infinite)and(hi.edge=hout)) or 
	  ((hi.cardinality=infinite)and(lo.edge=lout))
	 )
      then Sr:=-1
      else
      if Sr=0 then
         begin
            AlignUp(hi.exp,hi.mantissa,lo.exp,-lo.mantissa,EDiff,M0,M1,Dummy);
            MDiff:=M0+M1;
	    D:=GetReal(EDiff,MDiff);
	    if D < DCurr.RHalve[PC]*HalveThreshold
	    then {already narrowed enough dont bother} 
	       Change:=true {otherwise possible to terminate early}
	    else begin
	       DCurr.RHalve[PC]:=D;
               NewOuter(1);