with Text_Io;
package body Random_Number is

    M3 : constant := 97;  
    Denom : constant := 16777216.0;
    Init_C : constant := 362436.0 / Denom;
    Cd : constant := 7654321.0 / Denom;
    Cm : constant := 16777213.0 / Denom;
    subtype Range_1 is Integer range 0 .. M1 - 1;
    subtype Range_2 is Integer range 0 .. M2 - 1;
    subtype Range_3 is Integer range 1 .. M3;
    -- I, J, K must be in the range 1 to 178, and not all 1 while L may be
    -- any integer ranging from 0 to 168. If nonnegative values are assigned
    -- to I, J, K, L outside the specified ranges, the generator will still
    -- perform satifactorily but not produce exactly the same bit patterns
    -- on different computers.
    I, J, K : Range_1;
    Ni, Nj : Integer;
    L : Range_2;
    C : Float;
    U : array (Range_2) of Float;

    -- this procedure initializes the table for
    -- F (M3, Range_3'last / 3 + 1, - mod 1.0), a lagged-Fibonacci
    -- sequence generator, and produces values for the arithmetic sequence.
    procedure Start (New_I : Seed_Range_1 := Default_I;
                     New_J : Seed_Range_1 := Default_J;
                     New_K : Seed_Range_1 := Default_K;
                     New_L : Seed_Range_1 := Default_L) is  
        S, T : Float;
        M : Range_1;
    begin
        I := New_I;
        J := New_J;
        K := New_K;
        L := New_L;
        Ni := Range_3'Last;
        Nj := (Range_3'Last / 3) + 1;
        C := Init_C;
        for Ii in Range_3 loop
            S := 0.0;
            T := 0.5;
            for Jj in 1 .. 24 loop
                M := (((J * I) mod M1) * K) mod M1;
                I := J;
                J := K;
                K := M;
                L := (53 * L + 1) mod M2;
                if ((L * M) mod 64) >= 32 then
                    S := S + T;
                end if;  
                T := 0.5 * T;
            end loop;  
            U (Ii) := S;
        end loop;
    end Start;

    function Next return Float is
        Temp : Float;
    begin
        Temp := U (Ni) - U (Nj);
        if Temp < 0.0 then
            Temp := Temp + 1.0;
        end if;
        U (Ni) := Temp;
        Ni := Ni - 1;
        if Ni = 0 then
            Ni := Range_3'Last;
        end if;  
        Nj := Nj - 1;
        if Nj = 0 then
            Nj := Range_3'Last;
        end if;  
        C := C - Cd;
        if C < 0.0 then
            C := C + Cm;
        end if;
        Temp := Temp - C;
        if Temp < 0.0 then
            Temp := Temp + 1.0;
        end if;  
        return Temp;
    end Next;

begin
    Start;  -- intialize table u
end Random_Number;