## 2 Mathematical Basis of Montgomery Modular Multiplication over Binary Extension Field  
The binary extension field, GF(\(2^{k}\)), can be constructed using an irreducible polynomial \(F\) of the form \(F=\sum_{j=0}^{k}f_{j}\cdot\Gamma^{j}\), where \(0\leq j\leq k-1\) and \(f_{k}=f_{0}=1\). The addition of polynomials in GF(\(2^{k}\)) is straightforward and involves performing bitwise exclusive-OR (XOR) operations [20; 25]. However, the modular multiplication of polynomials within GF(\(2^{k}\)) presents additional complexities. This is because the intermediate results obtained during multiplication require a subsequent step of modular reduction using \(\Gamma^{k}=\sum_{j=0}^{k-1}f_{i}\cdot\Gamma^{j}\) to ensure that the resulting polynomial remains within GF(\(2^{k}\)) and has a degree lower than \(k\)[25]. By applying this reduction, the integrity of the binary extension field is preserved, enabling accurate and valid multiplication operations to be carried out effectively within GF(\(2^{k}\)). Let us consider \(\sigma\) and \(v\) as two elements in the binary extension field GF(\(2^{k}\)), and we aim to perform Montgomery Modular Multiplication (MMM) on them [20; 21]. We assume the existence of a special element \(\vartheta\) that satisfies the condition \(gcd(\vartheta,F)=1\). The Montgomery residues of \(\sigma\) and \(v\) are denoted by \(C\) and \(D\), respectively. The computation \(P=CD\vartheta^{-1}\bmod F=v\vartheta\bmod F=\sum_{j=0}^{k-1}p_{j}\cdot\Gamma ^{j}\) represents the MMM operation between the polynomials \(C=\sigma\vartheta\bmod F=\sum_{j=0}^{k-1}c_{j}\cdot\Gamma^{j}\) and \(D=v\ \vartheta\bmod F=\sum_{j=0}^{k-1}d_{j}\cdot\Gamma^{j}\)[20; 21]. To obtain the ultimate outcome, \(U\), the deployment of Montgomery multiplication with the inputs \(P\) and \(1\) becomes paramount. This transformative operation culminates in the equation \(U=P\vartheta^{-1}\bmod F=\sigma\ v\bmod F\). The utilization of Montgomery multiplication holds a distinct advantage, especially in a diverse range of applications characterized by repetitive multiplications. Notably, in the domains of elliptic curve point multiplication, inversion, and exponentiation, this technique assumes a central role. Its ability to seamlessly handle pre- and post-transformations sets it apart, rendering it a powerful tool for achieving efficient and accurate results. By selecting \(\vartheta=F^{\,(k-1)/2}\), the multiplication \(P=CD\vartheta^{-1}\bmod F\) can consistently be expressed in this form [21]. \[P=C(d_{0}+d_{1}F+\dots+d_{(k-1)/2}F^{\,(k-1)/2}+\dots+d_{k-1}F^{\,k-1})F^{\,-(k -1)/2}\bmod F \tag{1}\]  
Equation (1) can be alternatively presented in the following manner [21]:  
\[ \begin{split} P=& C(d_{(k-1)/2}+d_{(k+1)/2}F^{\,1}+ \dots+d_{(k-1)}F^{\,(k-1)/2})\\ &+C(d_{0}F^{\,-(k-1)/2}+d_{1}F^{\,-(k-3)/2}+\\ &\dots+d_{(k-3)/2}^{-1})\bmod F\end{split}  \tag{2}\]  
It is worth highlighting that, in the realm of practical applications, the common practice is to utilize odd values for \(k\). Therefore, when constructing the multiplier, the focus will primarily be on catering to the requirements and considerations specific to odd \(k\) values. Equation (2) can be split into two separate polynomials, denoted by A and B, which can be expressed as follows. The sum of A and B should construct Equation (2) [21]. \[A=Cd_{(k-1)}F^{\,(k-1)/2}+Cd_{(k-2)}F^{\,(k-3)/2}+\dots+Cd_{(k+1)/2}F^{\,1}+ Cd_{(k-1)/2}\bmod F \tag{3}\]  
\[B=Cd_{0}F^{\,-(k-1)/2}+Cd_{1}F^{\,-(k-3)/2}+\dots+Cd_{(k-5)/2}F^{\,-2}+Cd_{(k -3)/2}F^{\,-1}\bmod F \tag{4}\]  
To derive the recursive forms of Equations (3) and (4), we can rearrange them into the following organized manner [21]:  
\[ \begin{split} A=&(\dots((Cd_{(k-1)})F\bmod F+Cd_{ (k-2)})F\bmod F+\dots+\\ Cd_{(k+1)/2})F\bmod F+Cd_{(k-1)/2}\end{split}  \tag{5}\]  
\[ \begin{split} B=&(\dots(((Cd_{0})F^{\,-1}\bmod F+ Cd_{1})F^{\,-1}\bmod F+\dots+\\ Cd_{(k-5)/2}F^{\,-1}\bmod F+Cd_{(k-3)/2})F^{\,-1}\bmod F \end{split}