gcc/ada/libgnat/s-expmod.adb

   1 ------------------------------------------------------------------------------
   2 --                                                                          --
   3 --                         GNAT RUN-TIME COMPONENTS                         --
   4 --                                                                          --
   5 --                       S Y S T E M . E X P _ M O D                        --
   6 --                                                                          --
   7 --                                 B o d y                                  --
   8 --                                                                          --
   9 --          Copyright (C) 1992-2023, Free Software Foundation, Inc.         --
  10 --                                                                          --
  11 -- GNAT is free software;  you can  redistribute it  and/or modify it under --
  12 -- terms of the  GNU General Public License as published  by the Free Soft- --
  13 -- ware  Foundation;  either version 3,  or (at your option) any later ver- --
  14 -- sion.  GNAT is distributed in the hope that it will be useful, but WITH- --
  15 -- OUT ANY WARRANTY;  without even the  implied warranty of MERCHANTABILITY --
  16 -- or FITNESS FOR A PARTICULAR PURPOSE.                                     --
  17 --                                                                          --
  18 -- As a special exception under Section 7 of GPL version 3, you are granted --
  19 -- additional permissions described in the GCC Runtime Library Exception,   --
  20 -- version 3.1, as published by the Free Software Foundation.               --
  21 --                                                                          --
  22 -- You should have received a copy of the GNU General Public License and    --
  23 -- a copy of the GCC Runtime Library Exception along with this program;     --
  24 -- see the files COPYING3 and COPYING.RUNTIME respectively.  If not, see    --
  25 -- <http://www.gnu.org/licenses/>.                                          --
  26 --                                                                          --
  27 -- GNAT was originally developed  by the GNAT team at  New York University. --
  28 -- Extensive contributions were provided by Ada Core Technologies Inc.      --
  29 --                                                                          --
  30 ------------------------------------------------------------------------------
  31
  32 --  Preconditions, postconditions, ghost code, loop invariants and assertions
  33 --  in this unit are meant for analysis only, not for run-time checking, as it
  34 --  would be too costly otherwise. This is enforced by setting the assertion
  35 --  policy to Ignore.
  36
  37 pragma Assertion_Policy (Pre            => Ignore,
  38                          Post           => Ignore,
  39                          Ghost          => Ignore,
  40                          Loop_Invariant => Ignore,
  41                          Assert         => Ignore);
  42
  43 with Ada.Numerics.Big_Numbers.Big_Integers_Ghost;
  44 use Ada.Numerics.Big_Numbers.Big_Integers_Ghost;
  45
  46 package body System.Exp_Mod
  47   with SPARK_Mode
  48 is
  49    use System.Unsigned_Types;
  50
  51    --  Local lemmas
  52
  53    procedure Lemma_Add_Mod (X, Y : Big_Natural; B : Big_Positive)
  54    with
  55      Ghost,
  56      Post => (X + Y) mod B = ((X mod B) + (Y mod B)) mod B;
  57
  58    procedure Lemma_Exp_Expand (A : Big_Integer; Exp : Natural)
  59    with
  60      Ghost,
  61      Post =>
  62        (if Exp rem 2 = 0 then
  63           A ** Exp = A ** (Exp / 2) * A ** (Exp / 2)
  64         else
  65           A ** Exp = A ** (Exp / 2) * A ** (Exp / 2) * A);
  66
  67    procedure Lemma_Exp_Mod (A : Big_Natural; Exp : Natural; B : Big_Positive)
  68    with
  69      Ghost,
  70      Subprogram_Variant => (Decreases => Exp),
  71      Post => ((A mod B) ** Exp) mod B = (A ** Exp) mod B;
  72
  73    procedure Lemma_Mod_Ident (A : Big_Natural; B : Big_Positive)
  74    with
  75      Ghost,
  76      Pre => A < B,
  77      Post => A mod B = A;
  78
  79    procedure Lemma_Mod_Mod (A : Big_Integer; B : Big_Positive)
  80    with
  81      Ghost,
  82      Post => A mod B mod B = A mod B;
  83
  84    procedure Lemma_Mult_Div (X : Big_Natural; Y : Big_Positive)
  85    with
  86      Ghost,
  87      Post => X * Y / Y = X;
  88
  89    procedure Lemma_Mult_Mod (X, Y : Big_Natural; B : Big_Positive)
  90    with
  91      Ghost,
  92      --  The following subprogram variant can be added as soon as supported
  93      --  Subprogram_Variant => (Decreases => Y),
  94      Post => (X * Y) mod B = ((X mod B) * (Y mod B)) mod B;
  95
  96    -----------------------------
  97    -- Local lemma null bodies --
  98    -----------------------------
  99
 100    procedure Lemma_Mod_Ident (A : Big_Natural; B : Big_Positive) is null;
 101    procedure Lemma_Mod_Mod (A : Big_Integer; B : Big_Positive) is null;
 102    procedure Lemma_Mult_Div (X : Big_Natural; Y : Big_Positive) is null;
 103
 104    -------------------
 105    -- Lemma_Add_Mod --
 106    -------------------
 107
 108    procedure Lemma_Add_Mod (X, Y : Big_Natural; B : Big_Positive) is
 109
 110       procedure Lemma_Euclidean_Mod (Q, F, R : Big_Natural) with
 111         Pre  => F /= 0,
 112         Post => (Q * F + R) mod F = R mod F,
 113         Subprogram_Variant => (Decreases => Q);
 114
 115       -------------------------
 116       -- Lemma_Euclidean_Mod --
 117       -------------------------
 118
 119       procedure Lemma_Euclidean_Mod (Q, F, R : Big_Natural) is
 120       begin
 121          if Q > 0 then
 122             Lemma_Euclidean_Mod (Q - 1, F, R);
 123          end if;
 124       end Lemma_Euclidean_Mod;
 125
 126       --  Local variables
 127
 128       Left  : constant Big_Natural := (X + Y) mod B;
 129       Right : constant Big_Natural := ((X mod B) + (Y mod B)) mod B;
 130       XQuot : constant Big_Natural := X / B;
 131       YQuot : constant Big_Natural := Y / B;
 132       AQuot : constant Big_Natural := (X mod B + Y mod B) / B;
 133    begin
 134       if Y /= 0 and B > 1 then
 135          pragma Assert (X = XQuot * B + X mod B);
 136          pragma Assert (Y = YQuot * B + Y mod B);
 137          pragma Assert
 138            (Left = ((XQuot + YQuot) * B + X mod B + Y mod B) mod B);
 139          pragma Assert (X mod B + Y mod B = AQuot * B + Right);
 140          pragma Assert (Left = ((XQuot + YQuot + AQuot) * B + Right) mod B);
 141          Lemma_Euclidean_Mod (XQuot + YQuot + AQuot, B, Right);
 142          pragma Assert (Left = (Right mod B));
 143          pragma Assert (Left = Right);
 144       end if;
 145    end Lemma_Add_Mod;
 146
 147    ----------------------
 148    -- Lemma_Exp_Expand --
 149    ----------------------
 150
 151    procedure Lemma_Exp_Expand (A : Big_Integer; Exp : Natural) is
 152    begin
 153       if Exp rem 2 = 0 then
 154          pragma Assert (Exp = Exp / 2 + Exp / 2);
 155       else
 156          pragma Assert (Exp = Exp / 2 + Exp / 2 + 1);
 157          pragma Assert (A ** Exp = A ** (Exp / 2) * A ** (Exp / 2 + 1));
 158          pragma Assert (A ** (Exp / 2 + 1) = A ** (Exp / 2) * A);
 159          pragma Assert (A ** Exp = A ** (Exp / 2) * A ** (Exp / 2) * A);
 160       end if;
 161    end Lemma_Exp_Expand;
 162
 163    -------------------
 164    -- Lemma_Exp_Mod --
 165    -------------------
 166
 167    procedure Lemma_Exp_Mod (A : Big_Natural; Exp : Natural; B : Big_Positive)
 168    is
 169    begin
 170       if Exp /= 0 then
 171          declare
 172             Left  : constant Big_Integer := ((A mod B) ** Exp) mod B;
 173             Right : constant Big_Integer := (A ** Exp) mod B;
 174          begin
 175             Lemma_Mult_Mod (A mod B, (A mod B) ** (Exp - 1), B);
 176             Lemma_Mod_Mod (A, B);
 177             Lemma_Exp_Mod (A, Exp - 1, B);
 178             Lemma_Mult_Mod (A, A ** (Exp - 1), B);
 179             pragma Assert
 180               ((A mod B) * (A mod B) ** (Exp - 1) = (A mod B) ** Exp);
 181             pragma Assert (A * A ** (Exp - 1) = A ** Exp);
 182             pragma Assert (Left = Right);
 183          end;
 184       end if;
 185    end Lemma_Exp_Mod;
 186
 187    --------------------
 188    -- Lemma_Mult_Mod --
 189    --------------------
 190
 191    procedure Lemma_Mult_Mod (X, Y : Big_Natural; B : Big_Positive) is
 192       Left : constant Big_Natural := (X * Y) mod B;
 193       Right : constant Big_Natural := ((X mod B) * (Y mod B)) mod B;
 194    begin
 195       if Y /= 0 and B > 1 then
 196          Lemma_Add_Mod (X * (Y - 1), X, B);
 197          Lemma_Mult_Mod (X, Y - 1, B);
 198          Lemma_Mod_Mod (X, B);
 199          Lemma_Add_Mod ((X mod B) * ((Y - 1) mod B), X mod B, B);
 200          Lemma_Add_Mod (Y - 1, 1, B);
 201          pragma Assert (((Y - 1) mod B + 1) mod B = Y mod B);
 202          if (Y - 1) mod B + 1 < B then
 203             Lemma_Mod_Ident ((Y - 1) mod B + 1, B);
 204             Lemma_Mod_Mod ((X mod B) * (Y mod B), B);
 205             pragma Assert (Left = Right);
 206          else
 207             pragma Assert (Y mod B = 0);
 208             pragma Assert (Y / B * B = Y);
 209             pragma Assert ((X * Y) mod B = (X * Y) - (X * Y) / B * B);
 210             pragma Assert
 211               ((X * Y) mod B = (X * Y) - (X * (Y / B) * B) / B * B);
 212             Lemma_Mult_Div (X * (Y / B), B);
 213             pragma Assert (Left = 0);
 214             pragma Assert (Left = Right);
 215          end if;
 216       end if;
 217    end Lemma_Mult_Mod;
 218
 219    -----------------
 220    -- Exp_Modular --
 221    -----------------
 222
 223    function Exp_Modular
 224      (Left    : Unsigned;
 225       Modulus : Unsigned;
 226       Right   : Natural) return Unsigned
 227    is
 228       Result : Unsigned := 1;
 229       Factor : Unsigned := Left;
 230       Exp    : Natural := Right;
 231
 232       function Mult (X, Y : Unsigned) return Unsigned is
 233         (Unsigned (Long_Long_Unsigned (X) * Long_Long_Unsigned (Y)
 234                     mod Long_Long_Unsigned (Modulus)))
 235       with
 236         Pre => Modulus /= 0;
 237       --  Modular multiplication. Note that we can't take advantage of the
 238       --  compiler's circuit, because the modulus is not known statically.
 239
 240       --  Local ghost variables, functions and lemmas
 241
 242       M : constant Big_Positive := Big (Modulus) with Ghost;
 243
 244       function Equal_Modulo (X, Y : Big_Integer) return Boolean is
 245          (X mod M = Y mod M)
 246       with
 247         Ghost,
 248         Pre => Modulus /= 0;
 249
 250       procedure Lemma_Mult (X, Y : Unsigned)
 251       with
 252         Ghost,
 253         Post => Big (Mult (X, Y)) = (Big (X) * Big (Y)) mod M
 254           and then Big (Mult (X, Y)) < M;
 255
 256       procedure Lemma_Mult (X, Y : Unsigned) is null;
 257
 258       Rest : Big_Integer with Ghost;
 259       --  Ghost variable to hold Factor**Exp between Exp and Factor updates
 260
 261    begin
 262       pragma Assert (Modulus /= 1);
 263
 264       --  We use the standard logarithmic approach, Exp gets shifted right
 265       --  testing successive low order bits and Factor is the value of the
 266       --  base raised to the next power of 2.
 267
 268       --  Note: it is not worth special casing the cases of base values -1,0,+1
 269       --  since the expander does this when the base is a literal, and other
 270       --  cases will be extremely rare.
 271
 272       if Exp /= 0 then
 273          loop
 274             pragma Loop_Invariant (Exp > 0);
 275             pragma Loop_Invariant (Result < Modulus);
 276             pragma Loop_Invariant (Equal_Modulo
 277               (Big (Result) * Big (Factor) ** Exp, Big (Left) ** Right));
 278             pragma Loop_Variant (Decreases => Exp);
 279
 280             if Exp rem 2 /= 0 then
 281                pragma Assert
 282                  (Big (Factor) ** Exp
 283                   = Big (Factor) * Big (Factor) ** (Exp - 1));
 284                pragma Assert (Equal_Modulo
 285                  ((Big (Result) * Big (Factor)) * Big (Factor) ** (Exp - 1),
 286                   Big (Left) ** Right));
 287                pragma Assert (Big (Factor) >= 0);
 288                Lemma_Mult_Mod (Big (Result) * Big (Factor),
 289                                   Big (Factor) ** (Exp - 1),
 290                                   Big (Modulus));
 291                Lemma_Mult (Result, Factor);
 292
 293                Result := Mult (Result, Factor);
 294
 295                Lemma_Mod_Ident (Big (Result), Big (Modulus));
 296                Lemma_Mod_Mod (Big (Factor) ** (Exp - 1), Big (Modulus));
 297                Lemma_Mult_Mod (Big (Result),
 298                                   Big (Factor) ** (Exp - 1),
 299                                   Big (Modulus));
 300                pragma Assert (Equal_Modulo
 301                  (Big (Result) * Big (Factor) ** (Exp - 1),
 302                   Big (Left) ** Right));
 303                Lemma_Exp_Expand (Big (Factor), Exp - 1);
 304                pragma Assert (Exp / 2 = (Exp - 1) / 2);
 305             end if;
 306
 307             Lemma_Exp_Expand (Big (Factor), Exp);
 308
 309             Exp := Exp / 2;
 310             exit when Exp = 0;
 311
 312             Rest := Big (Factor) ** Exp;
 313             pragma Assert (Equal_Modulo
 314               (Big (Result) * (Rest * Rest), Big (Left) ** Right));
 315             Lemma_Exp_Mod (Big (Factor) * Big (Factor), Exp, Big (Modulus));
 316             pragma Assert
 317               ((Big (Factor) * Big (Factor)) ** Exp = Rest * Rest);
 318             pragma Assert (Equal_Modulo
 319               ((Big (Factor) * Big (Factor)) ** Exp,
 320                Rest * Rest));
 321             Lemma_Mult (Factor, Factor);
 322
 323             Factor := Mult (Factor, Factor);
 324
 325             Lemma_Mod_Mod (Rest * Rest, Big (Modulus));
 326             Lemma_Mod_Ident (Big (Result), Big (Modulus));
 327             Lemma_Mult_Mod (Big (Result), Rest * Rest, Big (Modulus));
 328             pragma Assert (Big (Factor) >= 0);
 329             Lemma_Mult_Mod (Big (Result), Big (Factor) ** Exp,
 330                                Big (Modulus));
 331             pragma Assert (Equal_Modulo
 332               (Big (Result) * Big (Factor) ** Exp, Big (Left) ** Right));
 333          end loop;
 334
 335          pragma Assert (Big (Result) = Big (Left) ** Right mod Big (Modulus));
 336       end if;
 337
 338       return Result;
 339
 340    end Exp_Modular;
 341
 342 end System.Exp_Mod;