-- file books\Ada\PIA10\website\perfect.txt

with Ada.Numerics.Big_Numbers.Big_Integers;
package Big renames Ada.Numerics.Big_Numbers.Big_Integers;

with Ada.Text_IO; use Ada.Text_IO;
with Big; use Big;
procedure Put_Num(N: in Big_Integer) is
begin
   Put(To_String(N, 0));
end Put_Num;

with Big; use Big;
with Ada.Text_IO, Ada.Integer_Text_IO;
use Ada.Text_IO, Ada.Integer_Text_IO;
with Put_Num;
procedure Perfection is

   Nop: constant := 30;

   Loop_Start, Loop_End: Integer;
   Level: Integer := 0;

   Is_Prime: Boolean;
   LL, MM: Big_Integer;

   Primes: constant array(1 .. Nop)  of Integer := 
   (3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,
    61,67,71,73,79,83,89,97,101,103,107,109,113,127);

procedure Lucas_Lehmer (Q: Integer;
                        Mersenne: out Big_Integer;
                        Lout: out Big_Integer;
                        Is_Prime: out Boolean) is

   M: constant Big_Integer := 2**Q - 1; 
   L: Big_Integer := 4;
begin
   for I in 3 .. Q loop
      if Level = 3 or I = Q then
         L := (L**2 - 2); -- no mod M for level 3 or if last iteration
      else
         L := (L**2 - 2) mod M;
      end if;
   end loop;
   Is_Prime := L mod M = 0;
   Lout := L;
   Mersenne := M;
end Lucas_Lehmer;

begin
   Put_Line("Welcome to the Mersenne Prime and Perfect Number program");
   loop
      Put_line("Give loop range in 1 to 30 for set of prime numbers");
      Put_Line("Start loop? ");  Get(Loop_Start);
      Put_Line("End loop? ");  Get(Loop_End);
      exit when Loop_Start in 1..Nop and Loop_End in 1..Nop and 
         Loop_End-Loop_Start >= 0;
      Put_Line("Sorry, not an acceptable range, please try again");
   end loop;
   loop
      Put_line("Level of detail required, answer 1, 2, or 3");
      Get(Level);
      exit when Level = 1 or Level = 2 or Level = 3;
   end loop;
   New_Line(2);
   Put_Line("Mersenne Primes");
   for I in Loop_Start .. Loop_End loop
      New_Line;
      Lucas_Lehmer(Primes(I), MM, LL, Is_Prime);
      Put(Primes(I), 4);  Put(" : ");
      Put_Num(MM);
      if Is_Prime then
         Put(" is prime");
         New_Line;
         Put_Num((MM+1)/2*MM);  Put(" is perfect");
         New_Line;
      else
         Put(" is not prime");
         New_Line;
      end if;
      New_Line;

      if Level /= 1 then
         Put("L is "); Put_Num(LL); New_Line(2); Put(" equals "); 
         Put_Num(MM); Put(" times "); New_Line(2); Put_Num(LL/MM);
         if not Is_Prime then
            New_Line; Put(" remainder = "); Put_Num(LL mod MM); 
            New_Line;
         end if;
       end if;
   end loop;
   Skip_Line(2);
end Perfection;