/* n+m is the size of p and q
330734
private static final int mula_small(byte[] p,byte[] q,int m,byte[] x,int n,int z) {
		int v=0;
		for (int i=0;i<n;++i) {
			v+=(q[i+m] & 0xFF)+z*(x[i] & 0xFF);
			p[i+m]=(byte)v;
			v>>=8;
		}
		return v;		
	}
415596
function ptsub(q, n) {
        var eps;

        if (n === 1 && q > 0.001 && q < 0.01) {
            eps = 1.0e-4;
        } else if (n === 2 && q < 0.0001) {
            eps = 1.0e-4;
        } else if (n === 1 && q < 0.001) {
            eps = 1.0e-2;
        } else {
            eps = 1.0e-5;
        }

        s = 10000.0;
        w = 0.0;
        while (true) {
            w = w + s;

            if (s <= eps) { return w; }
            if ((qe = 2.0 - p_t(n, w) * 2.0 - q) === 0.0) { return w; }
            if (qe < 0.0) {
                w = w - s;
                s = s / 10.0;
            }
        }
    }
1351315
function pnorm(qn) {
        var b = [1.570796288, 0.03706987906, -0.8364353589e-3,
                -0.2250947176e-3, 0.6841218299e-5, 0.5824238515e-5,
                -0.104527497e-5, 0.8360937017e-7, -0.3231081277e-8,
                0.3657763036e-10, 0.6936233982e-12],
            w1 = qn,
            w3 = -Math.log(4.0 * w1 * (1.0 - w1)),
            i = 1;

        if (qn < 0.0 || qn > 1.0) { return 0.0; }
        if (qn === 0.5) { return 0.0; }
        if (qn > 0.5) { w1 = 1.0 - w1; }

        w1 = b[0];
        for (i; i < 11; i++) {
            w1 += b[i] * Math.pow(w3, i);
        }

        if (qn > 0.5) { return Math.sqrt(w1 * w3); }

        return -Math.sqrt(w1 * w3);
    }
396901
function mpn_sub(r, n1, n2, l) {
	var c = 0;
	for (var i = 0; i < l; ++i) {
		r[i] = (c = (c >> 31) + n1[i] - n2[i]) & 0x7FFFFFFF;
	}
	return c >> 31;
}
348216
function mpn_unpack(n) {
	var k = n.length, words = new Uint32Array(Math.floor(k * 31 / 32)), l, lr = 0, w = 0, wr = 32;
	for (var i = 0, j = 0;;) {
		if (lr == 0) {
			if (i >= k) {
				if (w != 0) {
					words[words.length - ++j] = w >>> wr;
				}
				return words;
			}
			l = n[i++], lr = 31;
		}
		var s = Math.min(lr, wr);
		w = w >>> s | (l & (1 << s) - 1) << 32 - s, l >>>= s;
		if ((wr -= s) == 0) {
			words[words.length - ++j] = w;
			w = 0, wr = 32;
		}
		lr -= s;
	}
}
343052
def get_size(self, m):
        
        nrow, ncol = 0, 0
        if m[0] == 'F':
            nrow = self.n
        elif m[0] == 'G':
            nrow = self.m

        if m[1] == 'x':
            ncol = self.n
        elif m[1] == 'y':
            ncol = self.m

        return nrow, ncol
583306
function fp_sub(r, n1, n2, p, l) {
	if (mpn_sub(r, n1, n2, l)) {
		mpn_add(r, r, p, l);
	}
	return r;
}
116516
function mpn_dbl(r, n, l) {
	var c = 0;
	for (var i = 0; i < l; ++i) {
		r[i] = (c = c >>> 31 | n[i] << 1) & 0x7FFFFFFF;
	}
	return c >>> 31;
}
439279
def rmultinomial(n, p, size=None):
    
    # Leaving size=None as the default means return value is 1d array
    # if not specified-- nicer.

    # Single value for p:
    if len(np.shape(p)) == 1:
        return np.random.multinomial(n, p, size)

    # Multiple values for p:
    if np.isscalar(n):
        n = n * np.ones(np.shape(p)[0], dtype=np.int)
    out = np.empty(np.shape(p))
    for i in xrange(np.shape(p)[0]):
        out[i, :] = np.random.multinomial(n[i], p[i,:], size)
    return out
473171
public int[] add(int[] p, int[] q) {
    int len = Math.max(p.length, q.length);
    int[] result = new int[len];
    for (int i = 0; i < len; i++) {
      if (i < p.length && i < q.length) {
        result[i] = add(p[i], q[i]);
      } else if (i < p.length) {
        result[i] = p[i];
      } else {
        result[i] = q[i];
      }
    }
    return result;
  }
335074
void pass_n(int fi, 
	      float in[],  int in0, int istride,
	      float out[], int out0, int ostride,
	      int sign, int factor, int product){
    int i = 0, j = 0;
    int k, k1;

    int m = n / factor;
    int q = n / product;
    int p_1 = product / factor;
    int jump = (factor - 1) * p_1;

    int e, e1;

    for (i = 0; i < m; i++) {
      out[out0+ostride*i]   = in[in0+istride*i];
      out[out0+ostride*i+1] = in[in0+istride*i+1]; }

    for (e = 1; e < (factor - 1) / 2 + 1; e++) {
      for (i = 0; i < m; i++) {
	int idx = i + e * m;
	int idxc = i + (factor - e) * m;
	out[out0+ostride*idx]    = in[in0+istride*idx]   + in[in0+istride*idxc];
	out[out0+ostride*idx+1]  = in[in0+istride*idx+1] + in[in0+istride*idxc+1];
	out[out0+ostride*idxc]   = in[in0+istride*idx]   - in[in0+istride*idxc];
	out[out0+ostride*idxc+1] = in[in0+istride*idx+1] - in[in0+istride*idxc+1];  }}

    /* e = 0 */

    for (i=0; i<m; i++) {
      in[in0+istride*i]   = out[out0+ostride*i];
      in[in0+istride*i+1] = out[out0+ostride*i+1]; }

    for (e1 = 1; e1 < (factor - 1) / 2 + 1; e1++){
      for (i = 0; i < m; i++){
	in[in0+istride*i]   += out[out0+ostride*(i + e1*m)];
	in[in0+istride*i+1] += out[out0+ostride*(i + e1*m)+1]; }}

    float twiddl[] = twiddle[fi][q];

    for (e = 1; e < (factor-1)/2 + 1; e++){
      int idx = e;
      float w_real, w_imag;

      int em = e * m;
      int ecm = (factor - e) * m;

      for (i = 0; i < m; i++) {
	in[in0+istride*(i+em)]    = out[out0+ostride*i];
	in[in0+istride*(i+em)+1]  = out[out0+ostride*i+1];
	in[in0+istride*(i+ecm)]   = out[out0+ostride*i];
	in[in0+istride*(i+ecm)+1] = out[out0+ostride*i+1];
      }

      for (e1 = 1; e1 < (factor - 1) / 2 + 1; e1++) {
	if (idx == 0) {
	  w_real = 1;
	  w_imag = 0; }
	else {
	  w_real =       twiddl[2*(idx-1)];
	  w_imag = -sign*twiddl[2*(idx-1)+1]; }
	for (i = 0; i < m; i++) {
	  float ap = w_real * out[out0+ostride*(i + e1 * m)];
	  float am = w_imag * out[out0+ostride*(i + (factor - e1) *m)+1];

	  float bp = w_real * out[out0+ostride*(i + e1 * m)+1];
	  float bm = w_imag * out[out0+ostride*(i + (factor - e1) *m)];

	  in[in0+istride*(i + em)]    += (ap - am);
	  in[in0+istride*(i + em)+1]  += (bp + bm);
	  in[in0+istride*(i + ecm)]   += (ap + am);
	  in[in0+istride*(i + ecm)+1] += (bp - bm);
	}
	idx += e;
	idx %= factor;
      }
    }

    i = 0;
    j = 0;

    /* k = 0 */
    for (k1 = 0; k1 < p_1; k1++) {
      out[out0+ostride*k1]   = in[in0+istride*k1];
      out[out0+ostride*k1+1] = in[in0+istride*k1+1]; }

    for (e1 = 1; e1 < factor; e1++){
      for (k1 = 0; k1 < p_1; k1++){
	out[out0+ostride*(k1 + e1 * p_1)]   = in[in0+istride*(k1 + e1 * m)];
	out[out0+ostride*(k1 + e1 * p_1)+1] = in[in0+istride*(k1 + e1 * m)+1];
      }}

    i = p_1;
    j = product;

    for (k = 1; k < q; k++) {
      for (k1 = 0; k1 < p_1; k1++) {
	out[out0+ostride*j]   = in[in0+istride*i];
	out[out0+ostride*j+1] = in[in0+istride*i+1];
	i++;
	j++;
      }
      j += jump;
    }

    i = p_1;
    j = product;

    for (k = 1; k < q; k++) {
      twiddl = twiddle[fi][k];
      for (k1 = 0; k1 < p_1; k1++) {
	for (e1 = 1; e1 < factor; e1++) {
	  float x_real = in[in0+istride*(i + e1 * m)];
	  float x_imag = in[in0+istride*(i + e1 * m)+1];

	  float w_real =       twiddl[2*(e1-1)];
	  float w_imag = -sign*twiddl[2*(e1-1)+1];

	  out[out0+ostride*(j + e1 * p_1)]   = w_real * x_real - w_imag * x_imag;
	  out[out0+ostride*(j + e1 * p_1)+1] = w_real * x_imag + w_imag * x_real;
	}
 	i++;
	j++;
      }
      j += jump;
    }
  }
1110753
def AREA(a, b):
        """"""
        return cmp(b[0] * b[1], a[0] * a[1]) or cmp(b[1], a[1]) or cmp(b[0], a[0])
171717
function (m, n) {
      m = [m[0] >>> 16, m[0] & 0xffff, m[1] >>> 16, m[1] & 0xffff]
      n = [n[0] >>> 16, n[0] & 0xffff, n[1] >>> 16, n[1] & 0xffff]
      var o = [0, 0, 0, 0]
      o[3] += m[3] * n[3]
      o[2] += o[3] >>> 16
      o[3] &= 0xffff
      o[2] += m[2] * n[3]
      o[1] += o[2] >>> 16
      o[2] &= 0xffff
      o[2] += m[3] * n[2]
      o[1] += o[2] >>> 16
      o[2] &= 0xffff
      o[1] += m[1] * n[3]
      o[0] += o[1] >>> 16
      o[1] &= 0xffff
      o[1] += m[2] * n[2]
      o[0] += o[1] >>> 16
      o[1] &= 0xffff
      o[1] += m[3] * n[1]
      o[0] += o[1] >>> 16
      o[1] &= 0xffff
      o[0] += (m[0] * n[3]) + (m[1] * n[2]) + (m[2] * n[1]) + (m[3] * n[0])
      o[0] &= 0xffff
      return [(o[0] << 16) | o[1], (o[2] << 16) | o[3]]
    }
10617
function f1(d, m, r) {
    const t = m + d;
    const I = (t << r) | (t >>> (32 - r));
    return ((sBox[0][I >>> 24] ^ sBox[1][(I >>> 16) & 255]) - sBox[2][(I >>> 8) & 255]) + sBox[3][I & 255];
  }
1358570
def algorithm_u(ns, m):
    
    def visit(nv, av):
        ps = [[] for _ in range(m)]
        for jv in range(nv):
            ps[av[jv + 1]].append(ns[jv])
        return ps

    def f(mu, nu, sigma, n, a):
        if mu == 2:
            yield visit(n, a)
        else:
            for v in f(mu - 1, nu - 1, (mu + sigma) % 2, n, a):
                yield v
        if nu == mu + 1:
            a[mu] = mu - 1
            yield visit(n, a)
            while a[nu] > 0:
                a[nu] = a[nu] - 1
                yield visit(n, a)
        elif nu > mu + 1:
            if (mu + sigma) % 2 == 1:
                a[nu - 1] = mu - 1
            else:
                a[mu] = mu - 1
            if (a[nu] + sigma) % 2 == 1:
                for v in b(mu, nu - 1, 0, n, a):
                    yield v
            else:
                for v in f(mu, nu - 1, 0, n, a):
                    yield v
            while a[nu] > 0:
                a[nu] = a[nu] - 1
                if (a[nu] + sigma) % 2 == 1:
                    for v in b(mu, nu - 1, 0, n, a):
                        yield v
                else:
                    for v in f(mu, nu - 1, 0, n, a):
                        yield v

    def b(mu, nu, sigma, n, a):
        if nu == mu + 1:
            while a[nu] < mu - 1:
                yield visit(n, a)
                a[nu] = a[nu] + 1
            yield visit(n, a)
            a[mu] = 0
        elif nu > mu + 1:
            if (a[nu] + sigma) % 2 == 1:
                for v in f(mu, nu - 1, 0, n, a):
                    yield v
            else:
                for v in b(mu, nu - 1, 0, n, a):
                    yield v
            while a[nu] < mu - 1:
                a[nu] = a[nu] + 1
                if (a[nu] + sigma) % 2 == 1:
                    for v in f(mu, nu - 1, 0, n, a):
                        yield v
                else:
                    for v in b(mu, nu - 1, 0, n, a):
                        yield v
            if (mu + sigma) % 2 == 1:
                a[nu - 1] = 0
            else:
                a[mu] = 0
        if mu == 2:
            yield visit(n, a)
        else:
            for v in b(mu - 1, nu - 1, (mu + sigma) % 2, n, a):
                yield v

    ng = len(ns)
    ag = [0] * (ng + 1)
    for j in range(1, m + 1):
        ag[ng - m + j] = j - 1
    return f(m, ng, 0, ng, ag) if m > 1 else [[ns]]
339408
function(q,a,b,x,s,t){
        return __add(__rol(__add(__add(a,q),__add(x,t)),s),b);
    }
150248
function ecp_pubkey(Q, p, a, G, d, l) {
	var R = ecp_new(l);
	ecp_proj(Q, ecp_mul(R, d, G, a, p, l), p, l);
}