Modular (Binary) Exponentiation example

cbruce
Posts: 135
Joined: Sep 12, 2007 19:13
Location: Dallas, Texas

Modular (Binary) Exponentiation example

.
Modular Exponentiation is an operation that is used quite frequently
in areas such as cryptography. The Diffie-Hellman method is a prime
(pun intended) example.

I actually use the below function with the parameters shown for some
unique purposes... but it has limited use in most real world
programming; due to the fact that it can only handle the relatively
small values shown. Remember, something like the Diffie-Hellman
method is normally using 300 digit or larger numbers... that means
that DH calls for a bignum library.

So just consider this an educational example of how binary
exponentiation works.

Note:
If anybody wants to help me out by explaining why the intermediate
result inside the function is only updated when working with the ones
bits of the exponent - feel free! I wrote this a few years ago, and
every time I use it I wonder about that.

I've seen lots of examples on the web, but I've never seen code that
even comes close to commenting correctly about why the individual
steps are doing what they do... which is ridiculous considering that
there are so few steps involved and that it is a function that is in
great use. Hopefully, this will help someone else.

Code: Select all

`CONST DBGPRINT  = false' -------------------------------'    Modular Exponentiation' -------------------------------Function Modular_Exponentiation(base_num as ulong, exponent as ulongint, modulus as ulong) as ulong    ' Modular Exponentiation by Repeated Squaring - (Binary Exponentiation)    '    ' Modular Exponentiation - Wikipedia    ' https://en.wikipedia.org/wiki/Modular_exponentiation    '    ' We are able to directly calculate things, like a 32-bit number     ' raised to the power of a 64-bit number and then take the     ' modulo of that result by a 32-bit number, by utilizing a     ' numerical method called "repeated squaring" to break the     ' exponentiation down into a series of smaller exponentiations     ' that are manageable without the need for a BIGNUM library.    '        dim b1 as ulongint = base_num    dim e1 as ulongint = exponent    dim result as ulongint = 1'                                                ' This temp var is just so I can print which '                                                ' bit of the exponent I'm working on below.'                                                dim i_bit as Integer     = 0'                                                IF DBGPRINT THEN'                                                    print "INIT base             = "; b1'                                                    print "INIT exponent (priv)  = "; e1'                                                    print "INIT modulus          = "; modulus'                                                    print "INIT result           = "; result'                                                END IF    '    b1 = b1 MOD modulus'                                                IF DBGPRINT THEN'                                                    print "b1 MOD modulus        = "; b1'                                                    print'                                                END IF    '    ' Binary Exponentiation...     ' ... We're working with the exponent as a series of bits.     ' We will be doing a SHR on the exponent in every pass of the     ' loop in order to get at the exponent's next bit.    ' When all of the ones bits are gone (leaving exponent = 0),     ' then we are done with the exponentiation.    do while exponent > 0'                                                IF DBGPRINT THEN'                                                    print "DO exponent           = "; bin(exponent,32)'                                                END IF        ' Binary Exponentiation...         ' ... We're working with the exponent as a series of bits.         ' If bit zero of the remaining exponent bits is a one bit...        if (exponent and 1) then            ' I'm not a maths guy... I don't know why we are only             ' updating the result when working with the ones bits            ' of the exponent.            ' The equation just says to do this now.            result = (result * b1) MOD modulus'                                                IF DBGPRINT THEN'                                                    print "result*b1 MOD modulus = "; result'                                                END IF        end if'                                                IF DBGPRINT THEN'                                                    print "b1                    = "; b1'                                                END IF        ''                                                IF DBGPRINT THEN'                                                    ' We want to see both the square and the MOD results.'                                                    print "(b1 ^ 2)              = "; (b1 * b1)'                                                END IF        b1 = (b1 * b1) MOD modulus'                                                IF DBGPRINT THEN'                                                    ' We want to see both the square and the MOD results.'                                                    print "b1 MOD modulus        = "; b1'                                                END IF        '        ' Binary Exponentiation...         ' ... We're working with the exponent as a series of bits.         ' Shift the exponent right one bit to get at the next bit.        exponent = exponent shr 1'                                                IF DBGPRINT THEN'                                                    ' Track which bit of the exponent is next...'                                                    i_bit += 1'                                                    print i_bit; tab(3); " - exponent SHR     = "; bin(exponent,32)'                                                    print'                                                END IF    loop'                                                IF DBGPRINT THEN'                                                    print "RESULT ................ "; result; "  "; ucase(hex(result))'                                                    print'                                                    print'                                                END IF    '    return resultend function' ------------------------------------------------------------dim base_num as ulong    = 5dim exponent as ulongint = 2982377858398747' This is the largest 32-bit prime. It will result in the largest ' 32-bit result set.dim modulus as ulong     = 4294967291'dim result as ulong      = 0'result = Modular_Exponentiation(base_num, exponent, modulus)print "result = "; resultend`
dodicat
Posts: 5755
Joined: Jan 10, 2006 20:30
Location: Scotland

Re: Modular (Binary) Exponentiation example

Dunno cbruce, it just has to be an odd number.
I suppose the ones bits of the exponent is a fancy way of saying an odd number?
Took me a minute to get past that bit, and spent a little while with google, like yourself, couldn't find an easy answer.

Here is a string way.
Also solves the rosetta code question https://rosettacode.org/wiki/Modular_exponentiation#C without gmp or whatever.

Code: Select all

` Function _divide(n1 As String,n2 As String,decimal_places As Integer=10,dpflag As String="s") As String    Dim As String number=n1,divisor=n2    dpflag=Lcase(dpflag)    'For MOD    Dim As Integer modstop    If dpflag="mod" Then         If Len(n1)<Len(n2) Then Return n1        If Len(n1)=Len(n2) Then            If n1<n2 Then Return n1        End If        modstop=Len(n1)-Len(n2)+1    End If    If dpflag<>"mod" Then        If dpflag<>"s"  Then dpflag="raw"     End If    Dim runcount As Integer    '_______  LOOK UP TABLES ______________    Dim Qmod(0 To 19) As Ubyte    Dim bool(0 To 19) As Ubyte    For z As Integer=0 To 19        Qmod(z)=(z Mod 10+48)        bool(z)=(-(10>z))    Next z    Dim answer As String   'THE ANSWER STRING          '_______ SET THE DECIMAL WHERE IT SHOULD BE AT _______    Dim As String part1,part2    #macro set(decimal)    #macro insert(s,char,position)    If position > 0 And position <=Len(s) Then        part1=Mid(s,1,position-1)        part2=Mid(s,position)        s=part1+char+part2    End If    #endmacro    insert(answer,".",decpos)    answer=thepoint+zeros+answer    If dpflag="raw" Then        answer=Mid(answer,1,decimal_places)    End If    #endmacro    '______________________________________________    '__________ SPLIT A STRING ABOUT A CHARACTRR __________    Dim As String var1,var2    Dim pst As Integer    #macro split(stri,char,var1,var2)    pst=Instr(stri,char)    var1="":var2=""    If pst<>0 Then        var1=Rtrim(Mid(stri,1,pst),".")        var2=Ltrim(Mid(stri,pst),".")    Else        var1=stri    End If    #endmacro        #macro Removepoint(s)    split(s,".",var1,var2)    #endmacro    '__________ GET THE SIGN AND CLEAR THE -ve __________________    Dim sign As String    If Left(number,1)="-" Xor Left (divisor,1)="-" Then sign="-"    If Left(number,1)="-" Then  number=Ltrim(number,"-")    If Left (divisor,1)="-" Then divisor=Ltrim(divisor,"-")        'DETERMINE THE DECIMAL POSITION BEFORE THE DIVISION    Dim As Integer lennint,lenddec,lend,lenn,difflen    split(number,".",var1,var2)    lennint=Len(var1)    split(divisor,".",var1,var2)    lenddec=Len(var2)        If Instr(number,".") Then         Removepoint(number)        number=var1+var2    End If    If Instr(divisor,".") Then         Removepoint(divisor)        divisor=var1+var2    End If    Dim As Integer numzeros    numzeros=Len(number)    number=Ltrim(number,"0"):divisor=Ltrim (divisor,"0")    numzeros=numzeros-Len(number)    lend=Len(divisor):lenn=Len(number)    If lend>lenn Then difflen=lend-lenn    Dim decpos As Integer=lenddec+lennint-lend+2-numzeros 'THE POSITION INDICATOR    Dim _sgn As Byte=-Sgn(decpos)    If _sgn=0 Then _sgn=1    Dim As String thepoint=String(_sgn,".") 'DECIMAL AT START (IF)    Dim As String zeros=String(-decpos+1,"0")'ZEROS AT START (IF) e.g. .0009    If dpflag<>"mod" Then        If Len(zeros) =0 Then dpflag="s"    End If    Dim As Integer runlength    If Len(zeros) Then         runlength=decimal_places        answer=String(Len(zeros)+runlength+10,"0")        If dpflag="raw" Then             runlength=1            answer=String(Len(zeros)+runlength+10,"0")            If decimal_places>Len(zeros) Then                runlength=runlength+(decimal_places-Len(zeros))                answer=String(Len(zeros)+runlength+10,"0")            End If        End If            Else        decimal_places=decimal_places+decpos        runlength=decimal_places        answer=String(Len(zeros)+runlength+10,"0")    End If    '___________DECIMAL POSITION DETERMINED  _____________        'SET UP THE VARIABLES AND START UP CONDITIONS    number=number+String(difflen+decimal_places,"0")    Dim count As Integer    Dim temp As String    Dim copytemp As String    Dim topstring As String    Dim copytopstring As String    Dim As Integer lenf,lens    Dim As Ubyte takeaway,subtractcarry    Dim As Integer n3,diff    If Ltrim(divisor,"0")="" Then Return "Error :division by zero"       lens=Len(divisor)    topstring=Left(number,lend)    copytopstring=topstring    Do        count=0        Do            count=count+1            copytemp=temp                        Do                '___________________ QUICK SUBTRACTION loop _________________                                              lenf=Len(topstring)                If  lens<lenf=0 Then 'not                    If Lens>lenf Then                        temp= "done"                        Exit Do                    End If                    If divisor>topstring Then                         temp= "done"                        Exit Do                    End If                End If                                diff=lenf-lens                temp=topstring                subtractcarry=0                                For n3=lenf-1 To diff Step -1                    takeaway= topstring[n3]-divisor[n3-diff]+10-subtractcarry                    temp[n3]=Qmod(takeaway)                    subtractcarry=bool(takeaway)                Next n3                 If subtractcarry=0 Then Exit Do                If n3=-1 Then Exit Do                For n3=n3 To 0 Step -1                     takeaway= topstring[n3]-38-subtractcarry                    temp[n3]=Qmod(takeaway)                    subtractcarry=bool(takeaway)                    If subtractcarry=0 Then Exit Do                Next n3                Exit Do                            Loop 'single run            temp=Ltrim(temp,"0")            If temp="" Then temp= "0"            topstring=temp        Loop Until temp="done"        ' INDIVIDUAL CHARACTERS CARVED OFF ________________               runcount=runcount+1        If count=1 Then            topstring=copytopstring+Mid(number,lend+runcount,1)        Else            topstring=copytemp+Mid(number,lend+runcount,1)        End If        copytopstring=topstring        topstring=Ltrim(topstring,"0")        If dpflag="mod" Then            If runcount=modstop Then                 If topstring="" Then Return "0"                Return Mid(topstring,1,Len(topstring)-1)            End If        End If        answer[runcount-1]=count+47        If topstring="" And runcount>Len(n1)+1 Then            Exit Do        End If    Loop Until runcount=runlength+1        ' END OF RUN TO REQUIRED DECIMAL PLACES    set(decimal) 'PUT IN THE DECIMAL POINT    'THERE IS ALWAYS A DECIMAL POINT SOMEWHERE IN THE ANSWER    'NOW GET RID OF IT IF IT IS REDUNDANT    answer=Rtrim(answer,"0")    answer=Rtrim(answer,".")    answer=Ltrim(answer,"0")    If answer="" Then Return "0"    Return sign+answerEnd FunctionDim Shared As Integer _Mod(0 To 99),_Div(0 To 99)For z As Integer=0 To 99:_Mod(z)=(z Mod 10+48):_Div(z)=z\10:Next        Function Qmult(a As String,b As String) As String        Var flag=0,la = Len(a),lb = Len(b)        If Len(b)>Len(a) Then flag=1:Swap a,b:Swap la,lb        Dim As Ubyte n,carry,ai        Var c =String(la+lb,"0")        For i As Integer =la-1 To 0 Step -1            carry=0:ai=a[i]-48            For j As Integer =lb-1 To 0 Step -1                n = ai * (b[j]-48) + (c[i+j+1]-48) + carry                carry =_Div(n):c[i+j+1]=_Mod(n)            Next j            c[i]+=carry        Next i        If flag Then Swap a,b        Return Ltrim(c,"0")    End Function    '=======================================================================    #define mod_(a,b) _divide((a),(b),,"mod")    #define div_(a,b) iif(len((a))<len((b)),"0",_divide((a),(b),-2))        Function Modular_Exponentiation(base_num As String, exponent As String, modulus As String) As String        Dim b1 As String = base_num        Dim e1 As String = exponent        Dim As String result="1"        b1 = mod_(b1,modulus)        Do While e1 <> "0"            Var L=Len(e1)-1            If e1[L] And 1 Then                result=mod_(Qmult(result,b1),modulus)            End If            b1=mod_(qmult(b1,b1),modulus)            e1=div_(e1,"2")        Loop        Return result    End Function        Dim base_num As String    = "5"    Dim exponent As String = "2982377858398747"    Dim modulus As String     = "4294967291"    Var ans=Modular_Exponentiation(base_num,exponent,modulus)    Print    Print ans        base_num="2988348162058574136915891421498819466320163312926952423791023078876139"    exponent="2351399303373464486466122544523690094744975233415544072992656881240319"    modulus="10000000000000000000000000000000000000000"    ans=Modular_Exponentiation(base_num,exponent,modulus)    Print ans    Print "Done"    Sleep      `
cbruce
Posts: 135
Joined: Sep 12, 2007 19:13
Location: Dallas, Texas

Re: Modular (Binary) Exponentiation example

dodicat, that string version is a really cool one! Thanks!
integer
Posts: 378
Joined: Feb 01, 2007 16:54
Location: usa

Re: Modular (Binary) Exponentiation example

This is not a full (reasonable) explanation, but it is a start.

positional notation:
decimal: 70506 --> 7*10^4 + 0*10^3 + 5*10^2 + 0*10^1 + 6*10^0
7*10000 + 0*1000 + 8*100 + 0*10 + 6*1 == 70506

binary: 10101 --> 1*2^4 + 0*2^3 + 1*2^2 + 0*2^1 + 1*2^0
16 + 0 + 4 + 0 + 1 == 21

This is the SUM of the values.

Exponentiation is the PRODUCT of the positional values

21 decimal == 10101 binary

in decimal: 3^(21)
in binary: 3^(10101)

1 * 3^16 = 43046721
0 * 3^8 = 0 not used
1 * 3^4 = 81
0 * 3^2 = 0 not used
1 * 3^1 = 3

43046721 * 81 * 3 == 3^21 = 10460353203
cbruce
Posts: 135
Joined: Sep 12, 2007 19:13
Location: Dallas, Texas

Re: Modular (Binary) Exponentiation example

integer, I don't "grok" it, but at least I can see it and can start to wrap my head around it. Thanks!