src/big.rs - incubator-milagro-crypto-rust - Git at Google

 /*
 Licensed to the Apache Software Foundation (ASF) under one
 or more contributor license agreements.  See the NOTICE file
 distributed with this work for additional information
 regarding copyright ownership.  The ASF licenses this file
 to you under the Apache License, Version 2.0 (the
 "License"); you may not use this file except in compliance
 with the License.  You may obtain a copy of the License at

   http://www.apache.org/licenses/LICENSE-2.0

 Unless required by applicable law or agreed to in writing,
 software distributed under the License is distributed on an
 "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
 KIND, either express or implied.  See the License for the
 specific language governing permissions and limitations
 under the License.
 */

 use super::super::arch;
 use super::super::arch::Chunk;

 use super::super::arch::DChunk;

 use super::dbig::DBig;
 use rand::RAND;

 pub use super::rom::BASEBITS;
 pub use super::rom::MODBYTES;
 use std::cmp::Ordering;
 use std::fmt;

 pub const NLEN: usize = (1 + ((8 * MODBYTES - 1) / BASEBITS));
 pub const DNLEN: usize = 2 * NLEN;
 pub const BMASK: Chunk = ((1 << BASEBITS) - 1);
 pub const HBITS: usize = (BASEBITS / 2);
 pub const HMASK: Chunk = ((1 << HBITS) - 1);
 pub const NEXCESS: isize = (1 << ((arch::CHUNK) - BASEBITS - 1));
 pub const BIGBITS: usize = (MODBYTES * 8);

 #[derive(Copy)]
 pub struct Big {
     pub w: [Chunk; NLEN],
 }

 impl Clone for Big {
     fn clone(&self) -> Big {
         *self
     }
 }

 impl fmt::Display for Big {
     fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
         write!(f, "Big: [ {} ]", self.tostring())
     }
 }

 impl fmt::Debug for Big {
     fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
         write!(f, "Big: [ {} ]", self.tostring())
     }
 }

 impl PartialEq for Big {
     fn eq(&self, other: &Big) -> bool {
         Big::comp(self, other) == 0
     }
 }

 impl Ord for Big {
     fn cmp(&self, other: &Big) -> Ordering {
         let r = Big::comp(self, other);
         if r > 0 {
             return Ordering::Greater;
         }
         if r < 0 {
             return Ordering::Less;
         }
         Ordering::Equal
     }
 }

 impl Eq for Big {}

 impl PartialOrd for Big {
     fn partial_cmp(&self, other: &Big) -> Option<Ordering> {
         Some(self.cmp(other))
     }
 }

 impl Big {
     pub fn new() -> Big {
         Big { w: [0; NLEN] }
     }

     /// Convert a integer to a Big
     pub fn new_int(x: isize) -> Big {
         let mut s = Big::new();
         s.w[0] = x as Chunk;
         s
     }

     /// Takes an array of integers and converts to a Big
     pub fn new_ints(a: &[Chunk]) -> Big {
         let mut s = Big::new();
         for i in 0..NLEN {
             s.w[i] = a[i]
         }
         s
     }

     /// Creates a copy of the Big
     pub fn new_copy(y: &Big) -> Big {
         let mut s = Big::new();
         for i in 0..NLEN {
             s.w[i] = y.w[i]
         }
         s
     }

     pub fn new_dcopy(y: &DBig) -> Big {
         let mut s = Big::new();
         for i in 0..NLEN {
             s.w[i] = y.w[i]
         }
         s
     }

     pub fn get(&self, i: usize) -> Chunk {
         self.w[i]
     }

     pub fn set(&mut self, i: usize, x: Chunk) {
         self.w[i] = x;
     }

     pub fn xortop(&mut self, x: Chunk) {
         self.w[NLEN - 1] ^= x;
     }

     pub fn ortop(&mut self, x: Chunk) {
         self.w[NLEN - 1] |= x;
     }

     /// test for zero
     pub fn iszilch(&self) -> bool {
         for i in 0..NLEN {
             if self.w[i] != 0 {
                 return false;
             }
         }
         true
     }

     /// set to zero
     pub fn zero(&mut self) {
         for i in 0..NLEN {
             self.w[i] = 0
         }
     }

     /// Test for equal to one
     pub fn isunity(&self) -> bool {
         for i in 1..NLEN {
             if self.w[i] != 0 {
                 return false;
             }
         }
         if self.w[0] != 1 {
             return false;
         }
         true
     }

     /// Set to one
     pub fn one(&mut self) {
         self.w[0] = 1;
         for i in 1..NLEN {
             self.w[i] = 0;
         }
     }

     /// Copy from another Big
     pub fn copy(&mut self, x: &Big) {
         for i in 0..NLEN {
             self.w[i] = x.w[i]
         }
     }

     /// Copy from a DBig
     pub fn dcopy(&mut self, x: &DBig) {
         for i in 0..NLEN {
             self.w[i] = x.w[i]
         }
     }

     /// Get top and bottom half of =x*y+c+r
     pub fn muladd(a: Chunk, b: Chunk, c: Chunk, r: Chunk) -> (Chunk, Chunk) {
         let prod: DChunk = (a as DChunk) * (b as DChunk) + (c as DChunk) + (r as DChunk);
         let bot = (prod & (BMASK as DChunk)) as Chunk;
         let top = (prod >> BASEBITS) as Chunk;
         (top, bot)
     }

     /// normalise Big - force all digits < 2^BASEBITS
     pub fn norm(&mut self) -> Chunk {
         let mut carry = 0 as Chunk;
         for i in 0..NLEN - 1 {
             let d = self.w[i] + carry;
             self.w[i] = d & BMASK;
             carry = d >> BASEBITS;
         }
         self.w[NLEN - 1] += carry;
         (self.w[NLEN - 1] >> ((8 * MODBYTES) % BASEBITS)) as Chunk
     }

     /// Conditional swap of two bigs depending on d using XOR - no branches
     pub fn cswap(&mut self, b: &mut Big, d: isize) {
         let mut c = d as Chunk;
         c = !(c - 1);
         for i in 0..NLEN {
             let t = c & (self.w[i] ^ b.w[i]);
             self.w[i] ^= t;
             b.w[i] ^= t;
         }
     }

     pub fn cmove(&mut self, g: &Big, d: isize) {
         let b = -d as Chunk;
         for i in 0..NLEN {
             self.w[i] ^= (self.w[i] ^ g.w[i]) & b;
         }
     }

     /// Shift right by less than a word
     pub fn fshr(&mut self, k: usize) -> isize {
         let n = k;
         let w = self.w[0] & ((1 << n) - 1); // shifted out part
         for i in 0..NLEN - 1 {
             self.w[i] = (self.w[i] >> k) | ((self.w[i + 1] << (BASEBITS - n)) & BMASK);
         }
         self.w[NLEN - 1] = self.w[NLEN - 1] >> k;
         return w as isize;
     }

     /// general shift right
     pub fn shr(&mut self, k: usize) {
         let n = k % BASEBITS;
         let m = k / BASEBITS;
         for i in 0..NLEN - m - 1 {
             self.w[i] = (self.w[m + i] >> n) | ((self.w[m + i + 1] << (BASEBITS - n)) & BMASK)
         }
         self.w[NLEN - m - 1] = self.w[NLEN - 1] >> n;
         for i in NLEN - m..NLEN {
             self.w[i] = 0
         }
     }

     /// Shift right by less than a word
     pub fn fshl(&mut self, k: usize) -> isize {
         let n = k;
         self.w[NLEN - 1] = (self.w[NLEN - 1] << n) | (self.w[NLEN - 2] >> (BASEBITS - n));
         for i in (1..NLEN - 1).rev() {
             self.w[i] = ((self.w[i] << k) & BMASK) | (self.w[i - 1] >> (BASEBITS - n));
         }
         self.w[0] = (self.w[0] << n) & BMASK;
         // return excess - only used in ff.c
         (self.w[NLEN - 1] >> ((8 * MODBYTES) % BASEBITS)) as isize
     }

     /// General shift left
     pub fn shl(&mut self, k: usize) {
         let n = k % BASEBITS;
         let m = k / BASEBITS;

         self.w[NLEN - 1] = self.w[NLEN - 1 - m] << n;
         if NLEN >= m + 2 {
             self.w[NLEN - 1] |= self.w[NLEN - m - 2] >> (BASEBITS - n)
         }
         for i in (m + 1..NLEN - 1).rev() {
             self.w[i] = ((self.w[i - m] << n) & BMASK) | (self.w[i - m - 1] >> (BASEBITS - n));
         }
         self.w[m] = (self.w[0] << n) & BMASK;
         for i in 0..m {
             self.w[i] = 0
         }
     }

     /// Return number of bits
     pub fn nbits(&self) -> usize {
         let mut k = NLEN - 1;
         let mut s = Big::new_copy(&self);
         s.norm();
         while (k as isize) >= 0 && s.w[k] == 0 {
             k = k.wrapping_sub(1)
         }
         if (k as isize) < 0 {
             return 0;
         }
         let mut bts = BASEBITS * k;
         let mut c = s.w[k];
         while c != 0 {
             c /= 2;
             bts += 1;
         }
         bts
     }

     // Convert to Hex String
     pub fn tostring(&self) -> String {
         let mut s = String::new();
         let mut len = self.nbits();

         if len % 4 == 0 {
             len /= 4;
         } else {
             len /= 4;
             len += 1;
         }
         let mb = (MODBYTES * 2) as usize;
         if len < mb {
             len = mb
         }

         for i in (0..len).rev() {
             let mut b = Big::new_copy(&self);
             b.shr(i * 4);
             s = s + &format!("{:X}", b.w[0] & 15);
         }
         s
     }

     /// From Hex String
     pub fn fromstring(val: String) -> Big {
         let mut res = Big::new();
         let len = val.len();
         let op = &val[0..1];
         let n = u8::from_str_radix(op, 16).unwrap();
         res.w[0] += n as Chunk;
         for i in 1..len {
             res.shl(4);
             let op = &val[i..=i];
             let n = u8::from_str_radix(op, 16).unwrap();
             res.w[0] += n as Chunk;
         }
         res
     }

     // Self += r
     pub fn add(&mut self, r: &Big) {
         for i in 0..NLEN {
             self.w[i] += r.w[i]
         }
     }

     //? Bitwise OR
     pub fn or(&mut self, r: &Big) {
         for i in 0..NLEN {
             self.w[i] |= r.w[i]
         }
     }

     /// Self * 2
     pub fn dbl(&mut self) {
         for i in 0..NLEN {
             self.w[i] += self.w[i]
         }
     }

     /// Return self + x
     pub fn plus(&self, x: &Big) -> Big {
         let mut s = Big::new();
         for i in 0..NLEN {
             s.w[i] = self.w[i] + x.w[i];
         }
         s
     }

     /// Return self + 1
     pub fn inc(&mut self, x: isize) {
         self.norm();
         self.w[0] += x as Chunk;
     }

     ///  Return self - x
     pub fn minus(&self, x: &Big) -> Big {
         let mut d = Big::new();
         for i in 0..NLEN {
             d.w[i] = self.w[i] - x.w[i];
         }
         d
     }

     // self -= x
     pub fn sub(&mut self, x: &Big) {
         for i in 0..NLEN {
             self.w[i] -= x.w[i];
         }
     }

     // reverse subtract this=x-this
     pub fn rsub(&mut self, x: &Big) {
         for i in 0..NLEN {
             self.w[i] = x.w[i] - self.w[i]
         }
     }

     // self-=x, where x is int
     pub fn dec(&mut self, x: isize) {
         self.norm();
         self.w[0] -= x as Chunk;
     }

     // self*=x, where x is small int<NEXCESS
     pub fn imul(&mut self, c: isize) {
         for i in 0..NLEN {
             self.w[i] *= c as Chunk;
         }
     }

     // convert this Big to byte array
     pub fn tobytearray(&mut self, b: &mut [u8], n: usize) {
         let mut c = Big::new_copy(self);
         c.norm();

         for i in (0..(MODBYTES as usize)).rev() {
             b[i + n] = (c.w[0] & 0xff) as u8;
             c.fshr(8);
         }
     }

     // convert from byte array to Big
     pub fn frombytearray(b: &[u8], n: usize) -> Big {
         let mut m = Big::new();
         for i in 0..(MODBYTES as usize) {
             m.fshl(8);
             m.w[0] += (b[i + n] & 0xff) as Chunk;
         }
         m
     }

     pub fn tobytes(&mut self, b: &mut [u8]) {
         self.tobytearray(b, 0)
     }

     pub fn frombytes(b: &[u8]) -> Big {
         Big::frombytearray(b, 0)
     }

     // self*=x, where x is >NEXCESS
     pub fn pmul(&mut self, c: isize) -> Chunk {
         let mut carry = 0 as Chunk;
         for i in 0..NLEN {
             let ak = self.w[i];
             let tuple = Big::muladd(ak, c as Chunk, carry, 0 as Chunk);
             carry = tuple.0;
             self.w[i] = tuple.1;
         }
         carry
     }

     // self*=c and catch overflow in DBig
     pub fn pxmul(&mut self, c: isize) -> DBig {
         let mut m = DBig::new();
         let mut carry = 0 as Chunk;
         for j in 0..NLEN {
             let tuple = Big::muladd(self.w[j], c as Chunk, carry, m.w[j]);
             carry = tuple.0;
             m.w[j] = tuple.1;
         }
         m.w[NLEN] = carry;
         m
     }

     // divide by 3
     pub fn div3(&mut self) -> Chunk {
         let mut carry = 0 as Chunk;
         self.norm();
         let base = 1 << BASEBITS;
         for i in (0..NLEN).rev() {
             let ak = carry * base + self.w[i];
             self.w[i] = ak / 3;
             carry = ak % 3;
         }
         carry
     }

     // return a*b where result fits in a Big
     pub fn smul(a: &Big, b: &Big) -> Big {
         let mut c = Big::new();
         for i in 0..NLEN {
             let mut carry = 0 as Chunk;
             for j in 0..NLEN {
                 if i + j < NLEN {
                     let tuple = Big::muladd(a.w[i], b.w[j], carry, c.w[i + j]);
                     carry = tuple.0;
                     c.w[i + j] = tuple.1;
                 }
             }
         }
         c
     }

     // Compare a and b, return 0 if a==b, -1 if a<b, +1 if a>b. Inputs must be normalised
     pub fn comp(a: &Big, b: &Big) -> isize {
         for i in (0..NLEN).rev() {
             if a.w[i] == b.w[i] {
                 continue;
             }
             if a.w[i] > b.w[i] {
                 return 1;
             } else {
                 return -1;
             }
         }
         0
     }

     // set x = x mod 2^m
     pub fn mod2m(&mut self, m: usize) {
         let wd = m / BASEBITS;
         let bt = m % BASEBITS;
         let msk = (1 << bt) - 1;
         self.w[wd] &= msk;
         for i in wd + 1..NLEN {
             self.w[i] = 0
         }
     }

     // Arazi and Qi inversion mod 256
     pub fn invmod256(a: isize) -> isize {
         let mut t1: isize = 0;
         let mut c = (a >> 1) & 1;
         t1 += c;
         t1 &= 1;
         t1 = 2 - t1;
         t1 <<= 1;
         let mut u = t1 + 1;

         // i=2
         let mut b = a & 3;
         t1 = u * b;
         t1 >>= 2;
         c = (a >> 2) & 3;
         let mut t2 = (u * c) & 3;
         t1 += t2;
         t1 *= u;
         t1 &= 3;
         t1 = 4 - t1;
         t1 <<= 2;
         u += t1;

         // i=4
         b = a & 15;
         t1 = u * b;
         t1 >>= 4;
         c = (a >> 4) & 15;
         t2 = (u * c) & 15;
         t1 += t2;
         t1 *= u;
         t1 &= 15;
         t1 = 16 - t1;
         t1 <<= 4;
         u += t1;

         u
     }

     /// Return parity
     pub fn parity(&self) -> isize {
         (self.w[0] % 2) as isize
     }

     // Return n-th bit
     pub fn bit(&self, n: usize) -> isize {
         if (self.w[n / (BASEBITS as usize)] & (1 << (n % BASEBITS))) > 0 {
             return 1;
         }
         0
     }

     /// Return n last bits
     pub fn lastbits(&mut self, n: usize) -> isize {
         let msk = ((1 << n) - 1) as Chunk;
         self.norm();
         (self.w[0] & msk) as isize
     }

     /// a = 1/a mod 2^256. This is very fast!
     pub fn invmod2m(&mut self) {
         let mut u = Big::new();
         let mut b = Big::new();
         let mut c = Big::new();

         u.inc(Big::invmod256(self.lastbits(8)));

         let mut i = 8;
         while i < BIGBITS {
             u.norm();
             b.copy(self);
             b.mod2m(i);
             let mut t1 = Big::smul(&u, &b);
             t1.shr(i);
             c.copy(self);
             c.shr(i);
             c.mod2m(i);

             let mut t2 = Big::smul(&u, &c);
             t2.mod2m(i);
             t1.add(&t2);
             t1.norm();
             b = Big::smul(&t1, &u);
             t1.copy(&b);
             t1.mod2m(i);

             t2.one();
             t2.shl(i);
             t1.rsub(&t2);
             t1.norm();
             t1.shl(i);
             u.add(&t1);
             i <<= 1;
         }
         u.mod2m(BIGBITS);
         self.copy(&u);
         self.norm();
     }

     /// reduce self mod m
     pub fn rmod(&mut self, n: &Big) {
         let mut k = 0;
         let mut m = Big::new_copy(n);
         let mut r = Big::new();
         self.norm();
         if Big::comp(self, &m) < 0 {
             return;
         }
         loop {
             m.fshl(1);
             k += 1;
             if Big::comp(self, &m) < 0 {
                 break;
             }
         }

         while k > 0 {
             m.fshr(1);

             r.copy(self);
             r.sub(&m);
             r.norm();
             self.cmove(
                 &r,
                 (1 - ((r.w[NLEN - 1] >> (arch::CHUNK - 1)) & 1)) as isize,
             );
             k -= 1;
         }
     }

     /// Divide self by m
     pub fn div(&mut self, n: &Big) {
         let mut k = 0;
         self.norm();
         let mut e = Big::new_int(1);
         let mut b = Big::new_copy(self);
         let mut m = Big::new_copy(n);
         let mut r = Big::new();
         self.zero();

         while Big::comp(&b, &m) >= 0 {
             e.fshl(1);
             m.fshl(1);
             k += 1;
         }

         while k > 0 {
             m.fshr(1);
             e.fshr(1);

             r.copy(&b);
             r.sub(&m);
             r.norm();
             let d = (1 - ((r.w[NLEN - 1] >> (arch::CHUNK - 1)) & 1)) as isize;
             b.cmove(&r, d);
             r.copy(self);
             r.add(&e);
             r.norm();
             self.cmove(&r, d);
             k -= 1;
         }
     }

     /// Get 8*MODBYTES size random number
     pub fn random(rng: &mut RAND) -> Big {
         let mut m = Big::new();
         let mut j = 0;
         let mut r: u8 = 0;
         // generate random Big

         for _ in 0..8 * (MODBYTES as usize) {
             if j == 0 {
                 r = rng.getbyte()
             } else {
                 r >>= 1
             }

             let b = (r as Chunk) & 1;
             m.shl(1);
             m.w[0] += b;
             j += 1;
             j &= 7;
         }
         m
     }

     /// Create random Big in portable way, one bit at a time
     pub fn randomnum(q: &Big, rng: &mut RAND) -> Big {
         let mut d = DBig::new();
         let mut j = 0;
         let mut r: u8 = 0;
         let t = Big::new_copy(q);
         for _ in 0..2 * t.nbits() {
             if j == 0 {
                 r = rng.getbyte();
             } else {
                 r >>= 1
             }

             let b = (r as Chunk) & 1;
             d.shl(1);
             d.w[0] += b;
             j += 1;
             j &= 7;
         }
         let m = d.dmod(q);
         m
     }

     /// Jacobi Symbol (this/p). Returns 0, 1 or -1
     pub fn jacobi(&mut self, p: &Big) -> isize {
         let mut m: usize = 0;
         let mut t = Big::new();
         let mut x = Big::new();
         let mut n = Big::new();
         let zilch = Big::new();
         let one = Big::new_int(1);
         if p.parity() == 0 || Big::comp(self, &zilch) == 0 || Big::comp(p, &one) <= 0 {
             return 0;
         }
         self.norm();

         x.copy(self);
         n.copy(p);
         x.rmod(p);

         while Big::comp(&n, &one) > 0 {
             if Big::comp(&x, &zilch) == 0 {
                 return 0;
             }
             let n8 = n.lastbits(3) as usize;
             let mut k = 0;
             while x.parity() == 0 {
                 k += 1;
                 x.shr(1);
             }
             if k % 2 == 1 {
                 m += (n8 * n8 - 1) / 8
             }
             m += (n8 - 1) * ((x.lastbits(2) as usize) - 1) / 4;
             t.copy(&n);
             t.rmod(&x);
             n.copy(&x);
             x.copy(&t);
             m %= 2;
         }
         if m == 0 {
             return 1;
         }
         -1
     }

     /// self = 1/self mod p. Binary method
     pub fn invmodp(&mut self, p: &Big) {
         self.rmod(p);
         let mut u = Big::new_copy(self);
         let mut v = Big::new_copy(p);
         let mut x1 = Big::new_int(1);
         let mut x2 = Big::new();
         let mut t = Big::new();
         let one = Big::new_int(1);

         while (Big::comp(&u, &one) != 0) && (Big::comp(&v, &one) != 0) {
             while u.parity() == 0 {
                 u.fshr(1);
                 if x1.parity() != 0 {
                     x1.add(p);
                     x1.norm();
                 }
                 x1.fshr(1);
             }
             while v.parity() == 0 {
                 v.fshr(1);
                 if x2.parity() != 0 {
                     x2.add(p);
                     x2.norm();
                 }
                 x2.fshr(1);
             }
             if Big::comp(&u, &v) >= 0 {
                 u.sub(&v);
                 u.norm();
                 if Big::comp(&x1, &x2) >= 0 {
                     x1.sub(&x2)
                 } else {
                     t.copy(p);
                     t.sub(&x2);
                     x1.add(&t);
                 }
                 x1.norm();
             } else {
                 v.sub(&u);
                 v.norm();
                 if Big::comp(&x2, &x1) >= 0 {
                     x2.sub(&x1)
                 } else {
                     t.copy(p);
                     t.sub(&x1);
                     x2.add(&t);
                 }
                 x2.norm();
             }
         }
         if Big::comp(&u, &one) == 0 {
             self.copy(&x1)
         } else {
             self.copy(&x2)
         }
     }

     /// return a*b as DBig
     pub fn mul(a: &Big, b: &Big) -> DBig {
         let mut c = DBig::new();
         let rm = BMASK as DChunk;
         let rb = BASEBITS;

         let mut d: [DChunk; DNLEN] = [0; DNLEN];
         for i in 0..NLEN {
             d[i] = (a.w[i] as DChunk) * (b.w[i] as DChunk);
         }
         let mut s = d[0];
         let mut t = s;
         c.w[0] = (t & rm) as Chunk;
         let mut co = t >> rb;
         for k in 1..NLEN {
             s += d[k];
             t = co + s;
             for i in 1 + k / 2..=k {
                 t += ((a.w[i] - a.w[k - i]) as DChunk) * ((b.w[k - i] - b.w[i]) as DChunk)
             }
             c.w[k] = (t & rm) as Chunk;
             co = t >> rb;
         }
         for k in NLEN..2 * NLEN - 1 {
             s -= d[k - NLEN];
             t = co + s;
             let mut i = 1 + k / 2;
             while i < NLEN {
                 t += ((a.w[i] - a.w[k - i]) as DChunk) * ((b.w[k - i] - b.w[i]) as DChunk);
                 i += 1;
             }

             c.w[k] = (t & rm) as Chunk;
             co = t >> rb;
         }
         c.w[2 * NLEN - 1] = co as Chunk;
         c
     }

     /// return a^2 as DBig
     pub fn sqr(a: &Big) -> DBig {
         let mut c = DBig::new();
         let rm = BMASK as DChunk;
         let rb = BASEBITS;

         let mut t = (a.w[0] as DChunk) * (a.w[0] as DChunk);
         c.w[0] = (t & rm) as Chunk;
         let mut co = t >> rb;

         let mut j = 1;
         while j < NLEN - 1 {
             t = (a.w[j] as DChunk) * (a.w[0] as DChunk);
             for i in 1..(j + 1) / 2 {
                 t += (a.w[j - i] as DChunk) * (a.w[i] as DChunk);
             }
             t += t;
             t += co;
             c.w[j] = (t & rm) as Chunk;
             co = t >> rb;
             j += 1;
             t = (a.w[j] as DChunk) * (a.w[0] as DChunk);
             for i in 1..(j + 1) / 2 {
                 t += (a.w[j - i] as DChunk) * (a.w[i] as DChunk);
             }
             t += t;
             t += co;
             t += (a.w[j / 2] as DChunk) * (a.w[j / 2] as DChunk);
             c.w[j] = (t & rm) as Chunk;
             co = t >> rb;
             j += 1;
         }

         j = NLEN + (NLEN % 2) - 1;
         while j < DNLEN - 3 {
             t = (a.w[NLEN - 1] as DChunk) * (a.w[j + 1 - NLEN] as DChunk);
             for i in j + 2 - NLEN..(j + 1) / 2 {
                 t += (a.w[j - i] as DChunk) * (a.w[i] as DChunk);
             }
             t += t;
             t += co;
             c.w[j] = (t & rm) as Chunk;
             co = t >> rb;
             j += 1;
             t = (a.w[NLEN - 1] as DChunk) * (a.w[j + 1 - NLEN] as DChunk);
             for i in j + 2 - NLEN..(j + 1) / 2 {
                 t += (a.w[j - i] as DChunk) * (a.w[i] as DChunk);
             }
             t += t;
             t += co;
             t += (a.w[j / 2] as DChunk) * (a.w[j / 2] as DChunk);
             c.w[j] = (t & rm) as Chunk;
             co = t >> rb;
             j += 1;
         }

         t = (a.w[NLEN - 2] as DChunk) * (a.w[NLEN - 1] as DChunk);
         t += t;
         t += co;
         c.w[DNLEN - 3] = (t & rm) as Chunk;
         co = t >> rb;

         t = (a.w[NLEN - 1] as DChunk) * (a.w[NLEN - 1] as DChunk) + co;
         c.w[DNLEN - 2] = (t & rm) as Chunk;
         co = t >> rb;
         c.w[DNLEN - 1] = co as Chunk;

         c
     }

     /// Montegomery Reduction
     pub fn monty(md: &Big, mc: Chunk, d: &mut DBig) -> Big {
         let mut b = Big::new();
         let rm = BMASK as DChunk;
         let rb = BASEBITS;

         let mut dd: [DChunk; NLEN] = [0; NLEN];
         let mut v: [Chunk; NLEN] = [0; NLEN];

         b.zero();

         let mut t = d.w[0] as DChunk;
         v[0] = (((t & rm) as Chunk).wrapping_mul(mc)) & BMASK;
         t += (v[0] as DChunk) * (md.w[0] as DChunk);
         let mut c = (d.w[1] as DChunk) + (t >> rb);
         let mut s: DChunk = 0;
         for k in 1..NLEN {
             t = c + s + (v[0] as DChunk) * (md.w[k] as DChunk);
             let mut i = 1 + k / 2;
             while i < k {
                 t += ((v[k - i] - v[i]) as DChunk) * ((md.w[i] - md.w[k - i]) as DChunk);
                 i += 1;
             }
             v[k] = (((t & rm) as Chunk).wrapping_mul(mc)) & BMASK;
             t += (v[k] as DChunk) * (md.w[0] as DChunk);
             c = (d.w[k + 1] as DChunk) + (t >> rb);
             dd[k] = (v[k] as DChunk) * (md.w[k] as DChunk);
             s += dd[k];
         }

         for k in NLEN..2 * NLEN - 1 {
             t = c + s;
             let mut i = 1 + k / 2;
             while i < NLEN {
                 t += ((v[k - i] - v[i]) as DChunk) * ((md.w[i] - md.w[k - i]) as DChunk);
                 i += 1;
             }
             b.w[k - NLEN] = (t & rm) as Chunk;
             c = (d.w[k + 1] as DChunk) + (t >> rb);
             s -= dd[k + 1 - NLEN];
         }
         b.w[NLEN - 1] = (c & rm) as Chunk;
         b
     }

     pub fn ssn(r: &mut Big, a: &Big, m: &mut Big) -> isize {
         let n = NLEN - 1;
         m.w[0] = (m.w[0] >> 1) | ((m.w[1] << (BASEBITS - 1)) & BMASK);
         r.w[0] = a.w[0] - m.w[0];
         let mut carry = r.w[0] >> BASEBITS;
         r.w[0] &= BMASK;
         for i in 1..n {
             m.w[i] = (m.w[i] >> 1) | ((m.w[i + 1] << (BASEBITS - 1)) & BMASK);
             r.w[i] = a.w[i] - m.w[i] + carry;
             carry = r.w[i] >> BASEBITS;
             r.w[i] &= BMASK;
         }
         m.w[n] >>= 1;
         r.w[n] = a.w[n] - m.w[n] + carry;
         ((r.w[n] >> (arch::CHUNK - 1)) & 1) as isize
     }

     /// return a*b mod m
     pub fn modmul(a1: &Big, b1: &Big, m: &Big) -> Big {
         let mut a = Big::new_copy(a1);
         let mut b = Big::new_copy(b1);
         a.rmod(m);
         b.rmod(m);
         let mut d = Big::mul(&a, &b);
         d.dmod(m)
     }

     /// return a^2 mod m
     pub fn modsqr(a1: &Big, m: &Big) -> Big {
         let mut a = Big::new_copy(a1);
         a.rmod(m);
         let mut d = Big::sqr(&a);
         d.dmod(m)
     }

     /// return -a mod m
     pub fn modneg(a1: &Big, m: &Big) -> Big {
         let mut a = Big::new_copy(a1);
         a.rmod(m);
         m.minus(&a)
     }

     /// return this^e mod m
     pub fn powmod(&mut self, e1: &Big, m: &Big) -> Big {
         self.norm();
         let mut e = Big::new_copy(e1);
         e.norm();
         let mut a = Big::new_int(1);
         let mut z = Big::new_copy(&e);
         let mut s = Big::new_copy(self);
         loop {
             let bt = z.parity();
             z.fshr(1);
             if bt == 1 {
                 a = Big::modmul(&a, &s, m)
             }
             if z.iszilch() {
                 break;
             }
             s = Big::modsqr(&s, m);
         }
         a
     }
 }
	/*
	Licensed to the Apache Software Foundation (ASF) under one
	or more contributor license agreements. See the NOTICE file
	distributed with this work for additional information
	regarding copyright ownership. The ASF licenses this file
	to you under the Apache License, Version 2.0 (the
	"License"); you may not use this file except in compliance
	with the License. You may obtain a copy of the License at

	http://www.apache.org/licenses/LICENSE-2.0

	Unless required by applicable law or agreed to in writing,
	software distributed under the License is distributed on an
	"AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
	KIND, either express or implied. See the License for the
	specific language governing permissions and limitations
	under the License.
	*/

	use super::super::arch;
	use super::super::arch::Chunk;

	use super::super::arch::DChunk;

	use super::dbig::DBig;
	use rand::RAND;

	pub use super::rom::BASEBITS;
	pub use super::rom::MODBYTES;
	use std::cmp::Ordering;
	use std::fmt;

	pub const NLEN: usize = (1 + ((8 * MODBYTES - 1) / BASEBITS));
	pub const DNLEN: usize = 2 * NLEN;
	pub const BMASK: Chunk = ((1 << BASEBITS) - 1);
	pub const HBITS: usize = (BASEBITS / 2);
	pub const HMASK: Chunk = ((1 << HBITS) - 1);
	pub const NEXCESS: isize = (1 << ((arch::CHUNK) - BASEBITS - 1));
	pub const BIGBITS: usize = (MODBYTES * 8);

	#[derive(Copy)]
	pub struct Big {
	pub w: [Chunk; NLEN],
	}

	impl Clone for Big {
	fn clone(&self) -> Big {
	*self
	}
	}

	impl fmt::Display for Big {
	fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
	write!(f, "Big: [ {} ]", self.tostring())
	}
	}

	impl fmt::Debug for Big {
	fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
	write!(f, "Big: [ {} ]", self.tostring())
	}
	}

	impl PartialEq for Big {
	fn eq(&self, other: &Big) -> bool {
	Big::comp(self, other) == 0
	}
	}

	impl Ord for Big {
	fn cmp(&self, other: &Big) -> Ordering {
	let r = Big::comp(self, other);
	if r > 0 {
	return Ordering::Greater;
	}
	if r < 0 {
	return Ordering::Less;
	}
	Ordering::Equal
	}
	}

	impl Eq for Big {}

	impl PartialOrd for Big {
	fn partial_cmp(&self, other: &Big) -> Option<Ordering> {
	Some(self.cmp(other))
	}
	}

	impl Big {
	pub fn new() -> Big {
	Big { w: [0; NLEN] }
	}

	/// Convert a integer to a Big
	pub fn new_int(x: isize) -> Big {
	let mut s = Big::new();
	s.w[0] = x as Chunk;
	s
	}

	/// Takes an array of integers and converts to a Big
	pub fn new_ints(a: &[Chunk]) -> Big {
	let mut s = Big::new();
	for i in 0..NLEN {
	s.w[i] = a[i]
	}
	s
	}

	/// Creates a copy of the Big
	pub fn new_copy(y: &Big) -> Big {
	let mut s = Big::new();
	for i in 0..NLEN {
	s.w[i] = y.w[i]
	}
	s
	}

	pub fn new_dcopy(y: &DBig) -> Big {
	let mut s = Big::new();
	for i in 0..NLEN {
	s.w[i] = y.w[i]
	}
	s
	}

	pub fn get(&self, i: usize) -> Chunk {
	self.w[i]
	}

	pub fn set(&mut self, i: usize, x: Chunk) {
	self.w[i] = x;
	}

	pub fn xortop(&mut self, x: Chunk) {
	self.w[NLEN - 1] ^= x;
	}

	pub fn ortop(&mut self, x: Chunk) {
	self.w[NLEN - 1] \|= x;
	}

	/// test for zero
	pub fn iszilch(&self) -> bool {
	for i in 0..NLEN {
	if self.w[i] != 0 {
	return false;
	}
	}
	true
	}

	/// set to zero
	pub fn zero(&mut self) {
	for i in 0..NLEN {
	self.w[i] = 0
	}
	}

	/// Test for equal to one
	pub fn isunity(&self) -> bool {
	for i in 1..NLEN {
	if self.w[i] != 0 {
	return false;
	}
	}
	if self.w[0] != 1 {
	return false;
	}
	true
	}

	/// Set to one
	pub fn one(&mut self) {
	self.w[0] = 1;
	for i in 1..NLEN {
	self.w[i] = 0;
	}
	}

	/// Copy from another Big
	pub fn copy(&mut self, x: &Big) {
	for i in 0..NLEN {
	self.w[i] = x.w[i]
	}
	}

	/// Copy from a DBig
	pub fn dcopy(&mut self, x: &DBig) {
	for i in 0..NLEN {
	self.w[i] = x.w[i]
	}
	}

	/// Get top and bottom half of =x*y+c+r
	pub fn muladd(a: Chunk, b: Chunk, c: Chunk, r: Chunk) -> (Chunk, Chunk) {
	let prod: DChunk = (a as DChunk) * (b as DChunk) + (c as DChunk) + (r as DChunk);
	let bot = (prod & (BMASK as DChunk)) as Chunk;
	let top = (prod >> BASEBITS) as Chunk;
	(top, bot)
	}

	/// normalise Big - force all digits < 2^BASEBITS
	pub fn norm(&mut self) -> Chunk {
	let mut carry = 0 as Chunk;
	for i in 0..NLEN - 1 {
	let d = self.w[i] + carry;
	self.w[i] = d & BMASK;
	carry = d >> BASEBITS;
	}
	self.w[NLEN - 1] += carry;
	(self.w[NLEN - 1] >> ((8 * MODBYTES) % BASEBITS)) as Chunk
	}

	/// Conditional swap of two bigs depending on d using XOR - no branches
	pub fn cswap(&mut self, b: &mut Big, d: isize) {
	let mut c = d as Chunk;
	c = !(c - 1);
	for i in 0..NLEN {
	let t = c & (self.w[i] ^ b.w[i]);
	self.w[i] ^= t;
	b.w[i] ^= t;
	}
	}

	pub fn cmove(&mut self, g: &Big, d: isize) {
	let b = -d as Chunk;
	for i in 0..NLEN {
	self.w[i] ^= (self.w[i] ^ g.w[i]) & b;
	}
	}

	/// Shift right by less than a word
	pub fn fshr(&mut self, k: usize) -> isize {
	let n = k;
	let w = self.w[0] & ((1 << n) - 1); // shifted out part
	for i in 0..NLEN - 1 {
	self.w[i] = (self.w[i] >> k) \| ((self.w[i + 1] << (BASEBITS - n)) & BMASK);
	}
	self.w[NLEN - 1] = self.w[NLEN - 1] >> k;
	return w as isize;
	}

	/// general shift right
	pub fn shr(&mut self, k: usize) {
	let n = k % BASEBITS;
	let m = k / BASEBITS;
	for i in 0..NLEN - m - 1 {
	self.w[i] = (self.w[m + i] >> n) \| ((self.w[m + i + 1] << (BASEBITS - n)) & BMASK)
	}
	self.w[NLEN - m - 1] = self.w[NLEN - 1] >> n;
	for i in NLEN - m..NLEN {
	self.w[i] = 0
	}
	}

	/// Shift right by less than a word
	pub fn fshl(&mut self, k: usize) -> isize {
	let n = k;
	self.w[NLEN - 1] = (self.w[NLEN - 1] << n) \| (self.w[NLEN - 2] >> (BASEBITS - n));
	for i in (1..NLEN - 1).rev() {
	self.w[i] = ((self.w[i] << k) & BMASK) \| (self.w[i - 1] >> (BASEBITS - n));
	}
	self.w[0] = (self.w[0] << n) & BMASK;
	// return excess - only used in ff.c
	(self.w[NLEN - 1] >> ((8 * MODBYTES) % BASEBITS)) as isize
	}

	/// General shift left
	pub fn shl(&mut self, k: usize) {
	let n = k % BASEBITS;
	let m = k / BASEBITS;

	self.w[NLEN - 1] = self.w[NLEN - 1 - m] << n;
	if NLEN >= m + 2 {
	self.w[NLEN - 1] \|= self.w[NLEN - m - 2] >> (BASEBITS - n)
	}
	for i in (m + 1..NLEN - 1).rev() {
	self.w[i] = ((self.w[i - m] << n) & BMASK) \| (self.w[i - m - 1] >> (BASEBITS - n));
	}
	self.w[m] = (self.w[0] << n) & BMASK;
	for i in 0..m {
	self.w[i] = 0
	}
	}

	/// Return number of bits
	pub fn nbits(&self) -> usize {
	let mut k = NLEN - 1;
	let mut s = Big::new_copy(&self);
	s.norm();
	while (k as isize) >= 0 && s.w[k] == 0 {
	k = k.wrapping_sub(1)
	}
	if (k as isize) < 0 {
	return 0;
	}
	let mut bts = BASEBITS * k;
	let mut c = s.w[k];
	while c != 0 {
	c /= 2;
	bts += 1;
	}
	bts
	}

	// Convert to Hex String
	pub fn tostring(&self) -> String {
	let mut s = String::new();
	let mut len = self.nbits();

	if len % 4 == 0 {
	len /= 4;
	} else {
	len /= 4;
	len += 1;
	}
	let mb = (MODBYTES * 2) as usize;
	if len < mb {
	len = mb
	}

	for i in (0..len).rev() {
	let mut b = Big::new_copy(&self);
	b.shr(i * 4);
	s = s + &format!("{:X}", b.w[0] & 15);
	}
	s
	}

	/// From Hex String
	pub fn fromstring(val: String) -> Big {
	let mut res = Big::new();
	let len = val.len();
	let op = &val[0..1];
	let n = u8::from_str_radix(op, 16).unwrap();
	res.w[0] += n as Chunk;
	for i in 1..len {
	res.shl(4);
	let op = &val[i..=i];
	let n = u8::from_str_radix(op, 16).unwrap();
	res.w[0] += n as Chunk;
	}
	res
	}

	// Self += r
	pub fn add(&mut self, r: &Big) {
	for i in 0..NLEN {
	self.w[i] += r.w[i]
	}
	}

	//? Bitwise OR
	pub fn or(&mut self, r: &Big) {
	for i in 0..NLEN {
	self.w[i] \|= r.w[i]
	}
	}

	/// Self * 2
	pub fn dbl(&mut self) {
	for i in 0..NLEN {
	self.w[i] += self.w[i]
	}
	}

	/// Return self + x
	pub fn plus(&self, x: &Big) -> Big {
	let mut s = Big::new();
	for i in 0..NLEN {
	s.w[i] = self.w[i] + x.w[i];
	}
	s
	}

	/// Return self + 1
	pub fn inc(&mut self, x: isize) {
	self.norm();
	self.w[0] += x as Chunk;
	}

	/// Return self - x
	pub fn minus(&self, x: &Big) -> Big {
	let mut d = Big::new();
	for i in 0..NLEN {
	d.w[i] = self.w[i] - x.w[i];
	}
	d
	}

	// self -= x
	pub fn sub(&mut self, x: &Big) {
	for i in 0..NLEN {
	self.w[i] -= x.w[i];
	}
	}

	// reverse subtract this=x-this
	pub fn rsub(&mut self, x: &Big) {
	for i in 0..NLEN {
	self.w[i] = x.w[i] - self.w[i]
	}
	}

	// self-=x, where x is int
	pub fn dec(&mut self, x: isize) {
	self.norm();
	self.w[0] -= x as Chunk;
	}

	// self*=x, where x is small int<NEXCESS
	pub fn imul(&mut self, c: isize) {
	for i in 0..NLEN {
	self.w[i] *= c as Chunk;
	}
	}

	// convert this Big to byte array
	pub fn tobytearray(&mut self, b: &mut [u8], n: usize) {
	let mut c = Big::new_copy(self);
	c.norm();

	for i in (0..(MODBYTES as usize)).rev() {
	b[i + n] = (c.w[0] & 0xff) as u8;
	c.fshr(8);
	}
	}

	// convert from byte array to Big
	pub fn frombytearray(b: &[u8], n: usize) -> Big {
	let mut m = Big::new();
	for i in 0..(MODBYTES as usize) {
	m.fshl(8);
	m.w[0] += (b[i + n] & 0xff) as Chunk;
	}
	m
	}

	pub fn tobytes(&mut self, b: &mut [u8]) {
	self.tobytearray(b, 0)
	}

	pub fn frombytes(b: &[u8]) -> Big {
	Big::frombytearray(b, 0)
	}

	// self*=x, where x is >NEXCESS
	pub fn pmul(&mut self, c: isize) -> Chunk {
	let mut carry = 0 as Chunk;
	for i in 0..NLEN {
	let ak = self.w[i];
	let tuple = Big::muladd(ak, c as Chunk, carry, 0 as Chunk);
	carry = tuple.0;
	self.w[i] = tuple.1;
	}
	carry
	}

	// self*=c and catch overflow in DBig
	pub fn pxmul(&mut self, c: isize) -> DBig {
	let mut m = DBig::new();
	let mut carry = 0 as Chunk;
	for j in 0..NLEN {
	let tuple = Big::muladd(self.w[j], c as Chunk, carry, m.w[j]);
	carry = tuple.0;
	m.w[j] = tuple.1;
	}
	m.w[NLEN] = carry;
	m
	}

	// divide by 3
	pub fn div3(&mut self) -> Chunk {
	let mut carry = 0 as Chunk;
	self.norm();
	let base = 1 << BASEBITS;
	for i in (0..NLEN).rev() {
	let ak = carry * base + self.w[i];
	self.w[i] = ak / 3;
	carry = ak % 3;
	}
	carry
	}

	// return a*b where result fits in a Big
	pub fn smul(a: &Big, b: &Big) -> Big {
	let mut c = Big::new();
	for i in 0..NLEN {
	let mut carry = 0 as Chunk;
	for j in 0..NLEN {
	if i + j < NLEN {
	let tuple = Big::muladd(a.w[i], b.w[j], carry, c.w[i + j]);
	carry = tuple.0;
	c.w[i + j] = tuple.1;
	}
	}
	}
	c
	}

	// Compare a and b, return 0 if a==b, -1 if a<b, +1 if a>b. Inputs must be normalised
	pub fn comp(a: &Big, b: &Big) -> isize {
	for i in (0..NLEN).rev() {
	if a.w[i] == b.w[i] {
	continue;
	}
	if a.w[i] > b.w[i] {
	return 1;
	} else {
	return -1;
	}
	}
	0
	}

	// set x = x mod 2^m
	pub fn mod2m(&mut self, m: usize) {
	let wd = m / BASEBITS;
	let bt = m % BASEBITS;
	let msk = (1 << bt) - 1;
	self.w[wd] &= msk;
	for i in wd + 1..NLEN {
	self.w[i] = 0
	}
	}

	// Arazi and Qi inversion mod 256
	pub fn invmod256(a: isize) -> isize {
	let mut t1: isize = 0;
	let mut c = (a >> 1) & 1;
	t1 += c;
	t1 &= 1;
	t1 = 2 - t1;
	t1 <<= 1;
	let mut u = t1 + 1;

	// i=2
	let mut b = a & 3;
	t1 = u * b;
	t1 >>= 2;
	c = (a >> 2) & 3;
	let mut t2 = (u * c) & 3;
	t1 += t2;
	t1 *= u;
	t1 &= 3;
	t1 = 4 - t1;
	t1 <<= 2;
	u += t1;

	// i=4
	b = a & 15;
	t1 = u * b;
	t1 >>= 4;
	c = (a >> 4) & 15;
	t2 = (u * c) & 15;
	t1 += t2;
	t1 *= u;
	t1 &= 15;
	t1 = 16 - t1;
	t1 <<= 4;
	u += t1;

	u
	}

	/// Return parity
	pub fn parity(&self) -> isize {
	(self.w[0] % 2) as isize
	}

	// Return n-th bit
	pub fn bit(&self, n: usize) -> isize {
	if (self.w[n / (BASEBITS as usize)] & (1 << (n % BASEBITS))) > 0 {
	return 1;
	}
	0
	}

	/// Return n last bits
	pub fn lastbits(&mut self, n: usize) -> isize {
	let msk = ((1 << n) - 1) as Chunk;
	self.norm();
	(self.w[0] & msk) as isize
	}

	/// a = 1/a mod 2^256. This is very fast!
	pub fn invmod2m(&mut self) {
	let mut u = Big::new();
	let mut b = Big::new();
	let mut c = Big::new();

	u.inc(Big::invmod256(self.lastbits(8)));

	let mut i = 8;
	while i < BIGBITS {
	u.norm();
	b.copy(self);
	b.mod2m(i);
	let mut t1 = Big::smul(&u, &b);
	t1.shr(i);
	c.copy(self);
	c.shr(i);
	c.mod2m(i);

	let mut t2 = Big::smul(&u, &c);
	t2.mod2m(i);
	t1.add(&t2);
	t1.norm();
	b = Big::smul(&t1, &u);
	t1.copy(&b);
	t1.mod2m(i);

	t2.one();
	t2.shl(i);
	t1.rsub(&t2);
	t1.norm();
	t1.shl(i);
	u.add(&t1);
	i <<= 1;
	}
	u.mod2m(BIGBITS);
	self.copy(&u);
	self.norm();
	}

	/// reduce self mod m
	pub fn rmod(&mut self, n: &Big) {
	let mut k = 0;
	let mut m = Big::new_copy(n);
	let mut r = Big::new();
	self.norm();
	if Big::comp(self, &m) < 0 {
	return;
	}
	loop {
	m.fshl(1);
	k += 1;
	if Big::comp(self, &m) < 0 {
	break;
	}
	}

	while k > 0 {
	m.fshr(1);

	r.copy(self);
	r.sub(&m);
	r.norm();
	self.cmove(
	&r,
	(1 - ((r.w[NLEN - 1] >> (arch::CHUNK - 1)) & 1)) as isize,
	);
	k -= 1;
	}
	}

	/// Divide self by m
	pub fn div(&mut self, n: &Big) {
	let mut k = 0;
	self.norm();
	let mut e = Big::new_int(1);
	let mut b = Big::new_copy(self);
	let mut m = Big::new_copy(n);
	let mut r = Big::new();
	self.zero();

	while Big::comp(&b, &m) >= 0 {
	e.fshl(1);
	m.fshl(1);
	k += 1;
	}

	while k > 0 {
	m.fshr(1);
	e.fshr(1);

	r.copy(&b);
	r.sub(&m);
	r.norm();
	let d = (1 - ((r.w[NLEN - 1] >> (arch::CHUNK - 1)) & 1)) as isize;
	b.cmove(&r, d);
	r.copy(self);
	r.add(&e);
	r.norm();
	self.cmove(&r, d);
	k -= 1;
	}
	}

	/// Get 8*MODBYTES size random number
	pub fn random(rng: &mut RAND) -> Big {
	let mut m = Big::new();
	let mut j = 0;
	let mut r: u8 = 0;
	// generate random Big

	for _ in 0..8 * (MODBYTES as usize) {
	if j == 0 {
	r = rng.getbyte()
	} else {
	r >>= 1
	}

	let b = (r as Chunk) & 1;
	m.shl(1);
	m.w[0] += b;
	j += 1;
	j &= 7;
	}
	m
	}

	/// Create random Big in portable way, one bit at a time
	pub fn randomnum(q: &Big, rng: &mut RAND) -> Big {
	let mut d = DBig::new();
	let mut j = 0;
	let mut r: u8 = 0;
	let t = Big::new_copy(q);
	for _ in 0..2 * t.nbits() {
	if j == 0 {
	r = rng.getbyte();
	} else {
	r >>= 1
	}

	let b = (r as Chunk) & 1;
	d.shl(1);
	d.w[0] += b;
	j += 1;
	j &= 7;
	}
	let m = d.dmod(q);
	m
	}

	/// Jacobi Symbol (this/p). Returns 0, 1 or -1
	pub fn jacobi(&mut self, p: &Big) -> isize {
	let mut m: usize = 0;
	let mut t = Big::new();
	let mut x = Big::new();
	let mut n = Big::new();
	let zilch = Big::new();
	let one = Big::new_int(1);
	if p.parity() == 0 \|\| Big::comp(self, &zilch) == 0 \|\| Big::comp(p, &one) <= 0 {
	return 0;
	}
	self.norm();

	x.copy(self);
	n.copy(p);
	x.rmod(p);

	while Big::comp(&n, &one) > 0 {
	if Big::comp(&x, &zilch) == 0 {
	return 0;
	}
	let n8 = n.lastbits(3) as usize;
	let mut k = 0;
	while x.parity() == 0 {
	k += 1;
	x.shr(1);
	}
	if k % 2 == 1 {
	m += (n8 * n8 - 1) / 8
	}
	m += (n8 - 1) * ((x.lastbits(2) as usize) - 1) / 4;
	t.copy(&n);
	t.rmod(&x);
	n.copy(&x);
	x.copy(&t);
	m %= 2;
	}
	if m == 0 {
	return 1;
	}
	-1
	}

	/// self = 1/self mod p. Binary method
	pub fn invmodp(&mut self, p: &Big) {
	self.rmod(p);
	let mut u = Big::new_copy(self);
	let mut v = Big::new_copy(p);
	let mut x1 = Big::new_int(1);
	let mut x2 = Big::new();
	let mut t = Big::new();
	let one = Big::new_int(1);

	while (Big::comp(&u, &one) != 0) && (Big::comp(&v, &one) != 0) {
	while u.parity() == 0 {
	u.fshr(1);
	if x1.parity() != 0 {
	x1.add(p);
	x1.norm();
	}
	x1.fshr(1);
	}
	while v.parity() == 0 {
	v.fshr(1);
	if x2.parity() != 0 {
	x2.add(p);
	x2.norm();
	}
	x2.fshr(1);
	}
	if Big::comp(&u, &v) >= 0 {
	u.sub(&v);
	u.norm();
	if Big::comp(&x1, &x2) >= 0 {
	x1.sub(&x2)
	} else {
	t.copy(p);
	t.sub(&x2);
	x1.add(&t);
	}
	x1.norm();
	} else {
	v.sub(&u);
	v.norm();
	if Big::comp(&x2, &x1) >= 0 {
	x2.sub(&x1)
	} else {
	t.copy(p);
	t.sub(&x1);
	x2.add(&t);
	}
	x2.norm();
	}
	}
	if Big::comp(&u, &one) == 0 {
	self.copy(&x1)
	} else {
	self.copy(&x2)
	}
	}

	/// return a*b as DBig
	pub fn mul(a: &Big, b: &Big) -> DBig {
	let mut c = DBig::new();
	let rm = BMASK as DChunk;
	let rb = BASEBITS;

	let mut d: [DChunk; DNLEN] = [0; DNLEN];
	for i in 0..NLEN {
	d[i] = (a.w[i] as DChunk) * (b.w[i] as DChunk);
	}
	let mut s = d[0];
	let mut t = s;
	c.w[0] = (t & rm) as Chunk;
	let mut co = t >> rb;
	for k in 1..NLEN {
	s += d[k];
	t = co + s;
	for i in 1 + k / 2..=k {
	t += ((a.w[i] - a.w[k - i]) as DChunk) * ((b.w[k - i] - b.w[i]) as DChunk)
	}
	c.w[k] = (t & rm) as Chunk;
	co = t >> rb;
	}
	for k in NLEN..2 * NLEN - 1 {
	s -= d[k - NLEN];
	t = co + s;
	let mut i = 1 + k / 2;
	while i < NLEN {
	t += ((a.w[i] - a.w[k - i]) as DChunk) * ((b.w[k - i] - b.w[i]) as DChunk);
	i += 1;
	}

	c.w[k] = (t & rm) as Chunk;
	co = t >> rb;
	}
	c.w[2 * NLEN - 1] = co as Chunk;
	c
	}

	/// return a^2 as DBig
	pub fn sqr(a: &Big) -> DBig {
	let mut c = DBig::new();
	let rm = BMASK as DChunk;
	let rb = BASEBITS;

	let mut t = (a.w[0] as DChunk) * (a.w[0] as DChunk);
	c.w[0] = (t & rm) as Chunk;
	let mut co = t >> rb;

	let mut j = 1;
	while j < NLEN - 1 {
	t = (a.w[j] as DChunk) * (a.w[0] as DChunk);
	for i in 1..(j + 1) / 2 {
	t += (a.w[j - i] as DChunk) * (a.w[i] as DChunk);
	}
	t += t;
	t += co;
	c.w[j] = (t & rm) as Chunk;
	co = t >> rb;
	j += 1;
	t = (a.w[j] as DChunk) * (a.w[0] as DChunk);
	for i in 1..(j + 1) / 2 {
	t += (a.w[j - i] as DChunk) * (a.w[i] as DChunk);
	}
	t += t;
	t += co;
	t += (a.w[j / 2] as DChunk) * (a.w[j / 2] as DChunk);
	c.w[j] = (t & rm) as Chunk;
	co = t >> rb;
	j += 1;
	}

	j = NLEN + (NLEN % 2) - 1;
	while j < DNLEN - 3 {
	t = (a.w[NLEN - 1] as DChunk) * (a.w[j + 1 - NLEN] as DChunk);
	for i in j + 2 - NLEN..(j + 1) / 2 {
	t += (a.w[j - i] as DChunk) * (a.w[i] as DChunk);
	}
	t += t;
	t += co;
	c.w[j] = (t & rm) as Chunk;
	co = t >> rb;
	j += 1;
	t = (a.w[NLEN - 1] as DChunk) * (a.w[j + 1 - NLEN] as DChunk);
	for i in j + 2 - NLEN..(j + 1) / 2 {
	t += (a.w[j - i] as DChunk) * (a.w[i] as DChunk);
	}
	t += t;
	t += co;
	t += (a.w[j / 2] as DChunk) * (a.w[j / 2] as DChunk);
	c.w[j] = (t & rm) as Chunk;
	co = t >> rb;
	j += 1;
	}

	t = (a.w[NLEN - 2] as DChunk) * (a.w[NLEN - 1] as DChunk);
	t += t;
	t += co;
	c.w[DNLEN - 3] = (t & rm) as Chunk;
	co = t >> rb;

	t = (a.w[NLEN - 1] as DChunk) * (a.w[NLEN - 1] as DChunk) + co;
	c.w[DNLEN - 2] = (t & rm) as Chunk;
	co = t >> rb;
	c.w[DNLEN - 1] = co as Chunk;

	c
	}

	/// Montegomery Reduction
	pub fn monty(md: &Big, mc: Chunk, d: &mut DBig) -> Big {
	let mut b = Big::new();
	let rm = BMASK as DChunk;
	let rb = BASEBITS;

	let mut dd: [DChunk; NLEN] = [0; NLEN];
	let mut v: [Chunk; NLEN] = [0; NLEN];

	b.zero();

	let mut t = d.w[0] as DChunk;
	v[0] = (((t & rm) as Chunk).wrapping_mul(mc)) & BMASK;
	t += (v[0] as DChunk) * (md.w[0] as DChunk);
	let mut c = (d.w[1] as DChunk) + (t >> rb);
	let mut s: DChunk = 0;
	for k in 1..NLEN {
	t = c + s + (v[0] as DChunk) * (md.w[k] as DChunk);
	let mut i = 1 + k / 2;
	while i < k {
	t += ((v[k - i] - v[i]) as DChunk) * ((md.w[i] - md.w[k - i]) as DChunk);
	i += 1;
	}
	v[k] = (((t & rm) as Chunk).wrapping_mul(mc)) & BMASK;
	t += (v[k] as DChunk) * (md.w[0] as DChunk);
	c = (d.w[k + 1] as DChunk) + (t >> rb);
	dd[k] = (v[k] as DChunk) * (md.w[k] as DChunk);
	s += dd[k];
	}

	for k in NLEN..2 * NLEN - 1 {
	t = c + s;
	let mut i = 1 + k / 2;
	while i < NLEN {
	t += ((v[k - i] - v[i]) as DChunk) * ((md.w[i] - md.w[k - i]) as DChunk);
	i += 1;
	}
	b.w[k - NLEN] = (t & rm) as Chunk;
	c = (d.w[k + 1] as DChunk) + (t >> rb);
	s -= dd[k + 1 - NLEN];
	}
	b.w[NLEN - 1] = (c & rm) as Chunk;
	b
	}

	pub fn ssn(r: &mut Big, a: &Big, m: &mut Big) -> isize {
	let n = NLEN - 1;
	m.w[0] = (m.w[0] >> 1) \| ((m.w[1] << (BASEBITS - 1)) & BMASK);
	r.w[0] = a.w[0] - m.w[0];
	let mut carry = r.w[0] >> BASEBITS;
	r.w[0] &= BMASK;
	for i in 1..n {
	m.w[i] = (m.w[i] >> 1) \| ((m.w[i + 1] << (BASEBITS - 1)) & BMASK);
	r.w[i] = a.w[i] - m.w[i] + carry;
	carry = r.w[i] >> BASEBITS;
	r.w[i] &= BMASK;
	}
	m.w[n] >>= 1;
	r.w[n] = a.w[n] - m.w[n] + carry;
	((r.w[n] >> (arch::CHUNK - 1)) & 1) as isize
	}

	/// return a*b mod m
	pub fn modmul(a1: &Big, b1: &Big, m: &Big) -> Big {
	let mut a = Big::new_copy(a1);
	let mut b = Big::new_copy(b1);
	a.rmod(m);
	b.rmod(m);
	let mut d = Big::mul(&a, &b);
	d.dmod(m)
	}

	/// return a^2 mod m
	pub fn modsqr(a1: &Big, m: &Big) -> Big {
	let mut a = Big::new_copy(a1);
	a.rmod(m);
	let mut d = Big::sqr(&a);
	d.dmod(m)
	}

	/// return -a mod m
	pub fn modneg(a1: &Big, m: &Big) -> Big {
	let mut a = Big::new_copy(a1);
	a.rmod(m);
	m.minus(&a)
	}

	/// return this^e mod m
	pub fn powmod(&mut self, e1: &Big, m: &Big) -> Big {
	self.norm();
	let mut e = Big::new_copy(e1);
	e.norm();
	let mut a = Big::new_int(1);
	let mut z = Big::new_copy(&e);
	let mut s = Big::new_copy(self);
	loop {
	let bt = z.parity();
	z.fshr(1);
	if bt == 1 {
	a = Big::modmul(&a, &s, m)
	}
	if z.iszilch() {
	break;
	}
	s = Big::modsqr(&s, m);
	}
	a
	}
	}