DataMuseum.dk

Presents historical artifacts from the history of:

Rational R1000/400

This is an automatic "excavation" of a thematic subset of
artifacts from Datamuseum.dk's BitArchive.

See our Wiki for more about Rational R1000/400

Excavated with: AutoArchaeologist - Free & Open Source Software.


top - download

⟦53a092ec0⟧ Ada Source

    Length: 6144 (0x1800)
    Types: Ada Source
    Notes: 03_class, FILE, R1k_Segment, e3_tag, function Kf_Power, seg_0130d5, separate Generic_Elementary_Functions

Derivation

└─⟦8527c1e9b⟧ Bits:30000544 8mm tape, Rational 1000, Arrival backup of disks in PAM's R1000
    └─ ⟦cfc2e13cd⟧ »Space Info Vol 2« 
        └─⟦this⟧ 

E3 Source Code



separate (Generic_Elementary_Functions)

function Kf_Power (X, Y : Float_Type) return Float_Type is

-- On input, X and Y are floating-point value in Float_Type;
-- On output, the value of X**Y is returned.

   Result : Float_Type;

   Xx, Yy, R1, R2, Z1, Z2, W1, W2, Q, Y1, Y2, Tmp, F : Common_Float;
   N, M, J : Common_Int;
   L, Ll : Positive;

   Two_To : constant array (Common_Int range -3 .. 3) of Common_Float :=
      (0.125, 0.25, 0.5, 1.0, 2.0, 4.0, 8.0);
   Thirty_Two_By_Log2 : constant :=
      16#2E.2A8EC_A5705_FC2EE_FA1FF_B41A4_74FA2_3AD5E#;

   Large_Threshold : constant Common_Float :=
      4.0 * Common_Float (Float_Type'Machine_Emax) /
         Common_Float (Float_Type'Base'Epsilon);


begin

-- Filter out exceptional cases.

   Xx := Common_Float (X);
   Yy := Common_Float (Y);

   if (Xx = 0.0) then
      if (Y < 0.0) then
         raise Constraint_Error;
      elsif (Y > 0.0) then
         return (0.0);
      else
         raise Argument_Error;
      end if;
   end if;

   if (Xx < 0.0) then
      raise Argument_Error;
   end if;

   if (Xx = 1.0 or Y = 0.0) then
      return (1.0);
   end if;

   if (Yy = 1.0) then
      return (X);
   end if;

   if abs (Yy) >= Large_Threshold then
      if (Yy > 0.0) then
         -- Yy := Common_Float (Common_Float'Machine_Radix **
         --                     Common_Float'Machine_Emax);
         -- return (Float_Type (Yy * Yy * Yy));
         raise Constraint_Error; --pbk
      else
         return (Float_Type (0.0));
      end if;
   end if;

   L  := (2 * Common_Float'Machine_Mantissa) / 5;
   Ll := Common_Float'Machine_Mantissa;

-- Get log(X) to high accuracy

   Kp_Logpw (Xx, Z1, Z2);

   Y1 := Leading_Part (Yy, L);
   Y2 := Yy - Y1;

   W1 := Y1 * Z1;
   W2 := Y1 * Z2 + (Y2 * Z1 + Y2 * Z2);

   Z1 := Leading_Part (W1 + W2, Ll - 2);       -- force storage
   Z2 := W1 - Z1;
   Z2 := W2 + Z2;
-- Z1 + Z2 is Y*log(X) to high precision


-- Now get exp(Z1+Z2)

-- Perform argument reduction in Common_Float.

   N   := Common_Int (Z1 * Thirty_Two_By_Log2);
   Tmp := Common_Float (N) * Log2_By_32_Lead;
   if abs (Z1) >= abs (Tmp) then
      R1 := Z1 - Tmp;
   else
      Tmp := 0.5 * Tmp;
      R1  := (Z1 - Tmp) - Tmp;
   end if;
   R2 := Z2 - Common_Float (N) * Log2_By_32_Trail;
   J  := N mod 32;
   M  := (N - J) / 32;

-- Step 2. Approximation.

   Q := Kf_Em1sm (R1, R2);

-- Step 3. Function value reconstruction.

-- We now reconstruct the exponential of the input argument
-- so that Exp(Z1+Z2) = 2**M * (W1 + W2).
-- The order of the computation below must be strictly observed.

   F := Two_To_Jby32 (J, Lead) + Two_To_Jby32 (J, Trail);

   W1 := Two_To_Jby32 (J, Lead);
   W2 := Two_To_Jby32 (J, Trail) + F * Q;


   case Radix is

      when 2 =>
         Tmp := W1 + W2;

      when others =>
         J   := M rem 4;
         M   := (M - J) / 4;
         W1  := W1 * Two_To (J);
         W2  := W2 * Two_To (J);
         Tmp := W1 + W2;

   end case;

   Result := Float_Type (Scale (Tmp, M));

   return (Result);


end Kf_Power;


E3 Meta Data

    nblk1=5
    nid=4
    hdr6=8
        [0x00] rec0=25 rec1=00 rec2=01 rec3=010
        [0x01] rec0=2b rec1=00 rec2=02 rec3=00a
        [0x02] rec0=2a rec1=00 rec2=05 rec3=00c
        [0x03] rec0=10 rec1=00 rec2=03 rec3=000
        [0x04] rec0=00 rec1=00 rec2=00 rec3=008
    tail 0x2170e744682b151a9cf28 0x42a00066462061e03
Free Block Chain:
  0x4: 0000  00 00 01 0c 80 10 20 4a 20 20 20 3a 3d 20 4d 20  ┆       J   := M ┆