Finite Field Arithmetic

fz_lomul.adb


   1 ------------------------------------------------------------------------------
   2 ------------------------------------------------------------------------------
   3 -- This file is part of 'Finite Field Arithmetic', aka 'FFA'.               --
   4 --                                                                          --
   5 -- (C) 2019 Stanislav Datskovskiy ( www.loper-os.org )                      --
   6 -- http://wot.deedbot.org/17215D118B7239507FAFED98B98228A001ABFFC7.html     --
   7 --                                                                          --
   8 -- You do not have, nor can you ever acquire the right to use, copy or      --
   9 -- distribute this software ; Should you use this software for any purpose, --
  10 -- or copy and distribute it to anyone or in any manner, you are breaking   --
  11 -- the laws of whatever soi-disant jurisdiction, and you promise to         --
  12 -- continue doing so for the indefinite future. In any case, please         --
  13 -- always : read and understand any software ; verify any PGP signatures    --
  14 -- that you use - for any purpose.                                          --
  15 --                                                                          --
  16 -- See also http://trilema.com/2015/a-new-software-licensing-paradigm .     --
  17 ------------------------------------------------------------------------------
  18 ------------------------------------------------------------------------------
  19 
  20 with Words;    use Words;
  21 with Word_Ops; use Word_Ops;
  22 with W_Mul;    use W_Mul;
  23 with FZ_Arith; use FZ_Arith;
  24 with FZ_Mul;   use FZ_Mul;
  25 
  26 
  27 -- "Low Multiplication" computes only the bottom half of the product XY.
  28 -- Presently, it is used solely in Barrett's Modular Reduction.
  29 
  30 package body FZ_LoMul is
  31       
  32    -- Low-Only Comba's multiplier. (CAUTION: UNBUFFERED)
  33    procedure FZ_Low_Mul_Comba(X     : in  FZ;
  34                               Y     : in  FZ;
  35                               XY    : out FZ) is
  36       
  37       -- Words in each multiplicand (and also in the half-product)
  38       L : constant Word_Index := X'Length;
  39       
  40       -- 3-word Accumulator
  41       A2, A1, A0 : Word := 0;
  42       
  43    begin
  44       
  45       -- Compute the lower half of the Product, which is all we want:
  46       for N in 0 .. L - 1 loop
  47          
  48          -- Compute the Nth (indexed from zero) column of the Half-Product
  49          declare
  50             
  51             -- The outputs of a Word multiplication
  52             Lo, Hi : Word;
  53             
  54             -- Carry for the Accumulator addition
  55             C      : WBool;
  56             
  57             -- Sum for Accumulator addition
  58             Sum    : Word;
  59             
  60          begin
  61             
  62             -- For lower half of XY, will go from 0 to N
  63             -- For upper half of XY, will go from N - L + 1 to L - 1
  64             for j in 0 .. N loop
  65                
  66                -- Hi:Lo := j-th Word of X  *  (N - j)-th Word of Y
  67                Mul_Word(X(X'First + j),
  68                         Y(Y'First - j + N),
  69                         Lo, Hi);
  70                
  71                -- Now add Hi:Lo into the Accumulator:
  72                
  73                -- A0 += Lo; C := Carry
  74                Sum := A0 + Lo;
  75                C   := W_Carry(A0, Lo, Sum);
  76                A0  := Sum;
  77                
  78                -- A1 += Hi + C; C := Carry
  79                Sum := A1 + Hi + C;
  80                C   := W_Carry(A1, Hi, Sum);
  81                A1  := Sum;
  82                
  83                -- A2 += A2 + C
  84                A2  := A2 + C;
  85                
  86             end loop;
  87             
  88             -- We now have the Nth (indexed from zero) word of XY
  89             XY(XY'First + N) := A0;
  90             
  91             -- Right-Shift the Accumulator by one Word width:
  92             A0 := A1;
  93             A1 := A2;
  94             A2 := 0;
  95          end;
  96          
  97       end loop;
  98       
  99    end FZ_Low_Mul_Comba;
 100    
 101    
 102    -- Low-Only Multiplier. (CAUTION: UNBUFFERED)
 103    procedure Low_Mul(X  : in  FZ;
 104                      Y  : in  FZ;
 105                      XY : out FZ) is
 106       
 107       -- L is the wordness of a multiplicand. Guaranteed to be a power of two.
 108       L : constant Word_Count := X'Length;
 109       
 110       -- K is HALF of the length of a multiplicand.
 111       K : constant Word_Index := L / 2;
 112       
 113       -- A 'KSeg' is the same length as HALF of a multiplicand.
 114       subtype KSeg is FZ(1 .. K);
 115       
 116       -- The two K-sized variables of the half-product equation:
 117       LH, HL : KSeg;
 118       
 119       -- Bottom and Top K-sized halves of the multiplicand X.
 120       XLo        : KSeg  renames    X(  X'First       ..   X'Last - K );
 121       XHi        : KSeg  renames    X(  X'First + K   ..   X'Last     );
 122       
 123       -- Bottom and Top K-sized halves of the multiplicand Y.
 124       YLo        : KSeg  renames    Y(  Y'First       ..   Y'Last - K );
 125       YHi        : KSeg  renames    Y(  Y'First + K   ..   Y'Last     );
 126       
 127       -- Top K-sized half of the half-product XY.
 128       XYHi       : KSeg  renames   XY( XY'First + K   ..  XY'Last     );
 129       
 130       -- Carry from individual term additions.
 131       C          : WBool;
 132       pragma Unreferenced(C);
 133       
 134    begin
 135       
 136       -- Recurse to FULL-width multiplication: XY := XLo * YLo
 137       FZ_Multiply_Unbuffered(XLo, YLo, XY);
 138       
 139       -- Recurse to HALF-width multiplication: LH := XLo * YHi
 140       FZ_Low_Multiply_Unbuffered(XLo, YHi, LH);
 141       
 142       -- Recurse to HALF-width multiplication: HL := XHi * YLo
 143       FZ_Low_Multiply_Unbuffered(XHi, YLo, HL);
 144       
 145       -- XY += 2^(K * Bitness) * LH
 146       FZ_Add_D(X => XYHi, Y => LH, Overflow => C);
 147       
 148       -- XY += 2^(K * Bitness) * HL
 149       FZ_Add_D(X => XYHi, Y => HL, Overflow => C);
 150       
 151    end Low_Mul;
 152    -- CAUTION: Inlining prohibited for Low_Mul !
 153    
 154    
 155    -- Low-Only Multiplier. (CAUTION: UNBUFFERED)
 156    procedure FZ_Low_Multiply_Unbuffered(X     : in  FZ;
 157                                         Y     : in  FZ;
 158                                         XY    : out FZ) is
 159       
 160       -- The length of either multiplicand
 161       L : constant Word_Count := X'Length;
 162       
 163    begin
 164       
 165       if L <= Low_Mul_Thresh then
 166          
 167          -- Base case:
 168          FZ_Low_Mul_Comba(X, Y, XY);
 169          
 170       else
 171          
 172          -- Recursive case:
 173          Low_Mul(X, Y, XY);
 174          
 175       end if;
 176       
 177    end FZ_Low_Multiply_Unbuffered;
 178    
 179    
 180    -- Low-Only Multiplier. Preserves the inputs.
 181    procedure FZ_Low_Multiply_Buffered(X     : in  FZ;
 182                                       Y     : in  FZ;
 183                                       XY    : out FZ) is
 184       
 185       -- Product buffer.
 186       P : FZ(1 .. X'Length);
 187       
 188    begin
 189       
 190       FZ_Low_Multiply_Unbuffered(X, Y, P);
 191       
 192       XY := P;
 193       
 194    end FZ_Low_Multiply_Buffered;
 195    
 196 end FZ_LoMul;