diff --git a/jni/libzrtp/sources/bnlib/bn.h b/jni/libzrtp/sources/bnlib/bn.h
index 5cc80f0..f4dedfc 100644
--- a/jni/libzrtp/sources/bnlib/bn.h
+++ b/jni/libzrtp/sources/bnlib/bn.h
@@ -30,7 +30,6 @@
 	unsigned allocated;
 };
 
-#ifndef SWIG
 /*
  * User-supplied function: if non-NULL, this is called during long-running
  * computations.  You may put Yield() calls in here to give CPU time to
@@ -227,7 +226,6 @@
 	struct BnBasePrecomp const *pre1, struct BigNum const *exp1,
 	struct BnBasePrecomp const *pre2, struct BigNum const *exp2,
 	struct BigNum const *mod);
-#endif /* SWIF */
 
 #ifdef __cplusplus
 }
diff --git a/jni/libzrtp/sources/bnlib/bn32.c b/jni/libzrtp/sources/bnlib/bn32.c
index ee0d257..b0a3b85 100644
--- a/jni/libzrtp/sources/bnlib/bn32.c
+++ b/jni/libzrtp/sources/bnlib/bn32.c
@@ -13,7 +13,7 @@
 #define HAVE_CONFIG_H 0
 #endif
 #if HAVE_CONFIG_H
-#include <bnconfig.h>
+#include "bnconfig.h"
 #endif
 
 /*
diff --git a/jni/libzrtp/sources/bnlib/bn64.c b/jni/libzrtp/sources/bnlib/bn64.c
index 23cf185..6cb98ba 100644
--- a/jni/libzrtp/sources/bnlib/bn64.c
+++ b/jni/libzrtp/sources/bnlib/bn64.c
@@ -13,7 +13,7 @@
 #define HAVE_CONFIG_H 0
 #endif
 #if HAVE_CONFIG_H
-#include <bnconfig.h>
+#include "bnconfig.h"
 #endif
 
 /*
diff --git a/jni/libzrtp/sources/bnlib/bnprint.c b/jni/libzrtp/sources/bnlib/bnprint.c
index a407248..080e4a8 100644
--- a/jni/libzrtp/sources/bnlib/bnprint.c
+++ b/jni/libzrtp/sources/bnlib/bnprint.c
@@ -112,7 +112,12 @@
         getAsciiDigit(&d, radix, s[i]);
         bnMulQ(X, X, radix);
 
-        bnAddQ(X, d);
+        if(!neg ) {
+            bnAddQ(X, d);
+        }
+        else {
+            bnSubQ(X, d );
+        }
     }
-    return(neg);
+    return( 0 );
 }
diff --git a/jni/libzrtp/sources/bnlib/bnprint.h b/jni/libzrtp/sources/bnlib/bnprint.h
index b10393a..38830fb 100644
--- a/jni/libzrtp/sources/bnlib/bnprint.h
+++ b/jni/libzrtp/sources/bnlib/bnprint.h
@@ -10,13 +10,11 @@
 
 struct BigNum;
 
-#ifndef SWIG
 int bnPrint(FILE *f, char const *prefix, struct BigNum const *bn,
 	char const *suffix);
-#endif
 
 /**
- * Convert an ASCII string into a BigNum.
+ * Convert and ASCII string into a BigNum.
  *
  * This function converts an ASCII string into a Big number. If the first
  * character of the string is a minus sign the big number is a negative number.
diff --git a/jni/libzrtp/sources/bnlib/ec/curve25519-donna.c b/jni/libzrtp/sources/bnlib/ec/curve25519-donna.c
deleted file mode 100644
index de11280..0000000
--- a/jni/libzrtp/sources/bnlib/ec/curve25519-donna.c
+++ /dev/null
@@ -1,731 +0,0 @@
-/* Copyright 2008, Google Inc.
- * All rights reserved.
- *
- * Redistribution and use in source and binary forms, with or without
- * modification, are permitted provided that the following conditions are
- * met:
- *
- *     * Redistributions of source code must retain the above copyright
- * notice, this list of conditions and the following disclaimer.
- *     * Redistributions in binary form must reproduce the above
- * copyright notice, this list of conditions and the following disclaimer
- * in the documentation and/or other materials provided with the
- * distribution.
- *     * Neither the name of Google Inc. nor the names of its
- * contributors may be used to endorse or promote products derived from
- * this software without specific prior written permission.
- *
- * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
- * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
- * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
- * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
- * OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
- * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
- * LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
- * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
- * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
- * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
- * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
- *
- * curve25519-donna: Curve25519 elliptic curve, public key function
- *
- * http://code.google.com/p/curve25519-donna/
- *
- * Adam Langley <agl@imperialviolet.org>
- *
- * Derived from public domain C code by Daniel J. Bernstein <djb@cr.yp.to>
- *
- * More information about curve25519 can be found here
- *   http://cr.yp.to/ecdh.html
- *
- * djb's sample implementation of curve25519 is written in a special assembly
- * language called qhasm and uses the floating point registers.
- *
- * This is, almost, a clean room reimplementation from the curve25519 paper. It
- * uses many of the tricks described therein. Only the crecip function is taken
- * from the sample implementation.
- */
-
-#include <string.h>
-#include <stdint.h>
-
-#ifdef _MSC_VER
-#define inline __inline
-#endif
-
-typedef uint8_t u8;
-typedef int32_t s32;
-typedef int64_t limb;
-
-/* Field element representation:
- *
- * Field elements are written as an array of signed, 64-bit limbs, least
- * significant first. The value of the field element is:
- *   x[0] + 2^26·x[1] + x^51·x[2] + 2^102·x[3] + ...
- *
- * i.e. the limbs are 26, 25, 26, 25, ... bits wide.
- */
-
-/* Sum two numbers: output += in */
-static void fsum(limb *output, const limb *in) {
-  unsigned i;
-  for (i = 0; i < 10; i += 2) {
-    output[0+i] = (output[0+i] + in[0+i]);
-    output[1+i] = (output[1+i] + in[1+i]);
-  }
-}
-
-/* Find the difference of two numbers: output = in - output
- * (note the order of the arguments!)
- */
-static void fdifference(limb *output, const limb *in) {
-  unsigned i;
-  for (i = 0; i < 10; ++i) {
-    output[i] = (in[i] - output[i]);
-  }
-}
-
-/* Multiply a number by a scalar: output = in * scalar */
-static void fscalar_product(limb *output, const limb *in, const limb scalar) {
-  unsigned i;
-  for (i = 0; i < 10; ++i) {
-    output[i] = in[i] * scalar;
-  }
-}
-
-/* Multiply two numbers: output = in2 * in
- *
- * output must be distinct to both inputs. The inputs are reduced coefficient
- * form, the output is not.
- */
-static void fproduct(limb *output, const limb *in2, const limb *in) {
-  output[0] =       ((limb) ((s32) in2[0])) * ((s32) in[0]);
-  output[1] =       ((limb) ((s32) in2[0])) * ((s32) in[1]) +
-                    ((limb) ((s32) in2[1])) * ((s32) in[0]);
-  output[2] =  2 *  ((limb) ((s32) in2[1])) * ((s32) in[1]) +
-                    ((limb) ((s32) in2[0])) * ((s32) in[2]) +
-                    ((limb) ((s32) in2[2])) * ((s32) in[0]);
-  output[3] =       ((limb) ((s32) in2[1])) * ((s32) in[2]) +
-                    ((limb) ((s32) in2[2])) * ((s32) in[1]) +
-                    ((limb) ((s32) in2[0])) * ((s32) in[3]) +
-                    ((limb) ((s32) in2[3])) * ((s32) in[0]);
-  output[4] =       ((limb) ((s32) in2[2])) * ((s32) in[2]) +
-               2 * (((limb) ((s32) in2[1])) * ((s32) in[3]) +
-                    ((limb) ((s32) in2[3])) * ((s32) in[1])) +
-                    ((limb) ((s32) in2[0])) * ((s32) in[4]) +
-                    ((limb) ((s32) in2[4])) * ((s32) in[0]);
-  output[5] =       ((limb) ((s32) in2[2])) * ((s32) in[3]) +
-                    ((limb) ((s32) in2[3])) * ((s32) in[2]) +
-                    ((limb) ((s32) in2[1])) * ((s32) in[4]) +
-                    ((limb) ((s32) in2[4])) * ((s32) in[1]) +
-                    ((limb) ((s32) in2[0])) * ((s32) in[5]) +
-                    ((limb) ((s32) in2[5])) * ((s32) in[0]);
-  output[6] =  2 * (((limb) ((s32) in2[3])) * ((s32) in[3]) +
-                    ((limb) ((s32) in2[1])) * ((s32) in[5]) +
-                    ((limb) ((s32) in2[5])) * ((s32) in[1])) +
-                    ((limb) ((s32) in2[2])) * ((s32) in[4]) +
-                    ((limb) ((s32) in2[4])) * ((s32) in[2]) +
-                    ((limb) ((s32) in2[0])) * ((s32) in[6]) +
-                    ((limb) ((s32) in2[6])) * ((s32) in[0]);
-  output[7] =       ((limb) ((s32) in2[3])) * ((s32) in[4]) +
-                    ((limb) ((s32) in2[4])) * ((s32) in[3]) +
-                    ((limb) ((s32) in2[2])) * ((s32) in[5]) +
-                    ((limb) ((s32) in2[5])) * ((s32) in[2]) +
-                    ((limb) ((s32) in2[1])) * ((s32) in[6]) +
-                    ((limb) ((s32) in2[6])) * ((s32) in[1]) +
-                    ((limb) ((s32) in2[0])) * ((s32) in[7]) +
-                    ((limb) ((s32) in2[7])) * ((s32) in[0]);
-  output[8] =       ((limb) ((s32) in2[4])) * ((s32) in[4]) +
-               2 * (((limb) ((s32) in2[3])) * ((s32) in[5]) +
-                    ((limb) ((s32) in2[5])) * ((s32) in[3]) +
-                    ((limb) ((s32) in2[1])) * ((s32) in[7]) +
-                    ((limb) ((s32) in2[7])) * ((s32) in[1])) +
-                    ((limb) ((s32) in2[2])) * ((s32) in[6]) +
-                    ((limb) ((s32) in2[6])) * ((s32) in[2]) +
-                    ((limb) ((s32) in2[0])) * ((s32) in[8]) +
-                    ((limb) ((s32) in2[8])) * ((s32) in[0]);
-  output[9] =       ((limb) ((s32) in2[4])) * ((s32) in[5]) +
-                    ((limb) ((s32) in2[5])) * ((s32) in[4]) +
-                    ((limb) ((s32) in2[3])) * ((s32) in[6]) +
-                    ((limb) ((s32) in2[6])) * ((s32) in[3]) +
-                    ((limb) ((s32) in2[2])) * ((s32) in[7]) +
-                    ((limb) ((s32) in2[7])) * ((s32) in[2]) +
-                    ((limb) ((s32) in2[1])) * ((s32) in[8]) +
-                    ((limb) ((s32) in2[8])) * ((s32) in[1]) +
-                    ((limb) ((s32) in2[0])) * ((s32) in[9]) +
-                    ((limb) ((s32) in2[9])) * ((s32) in[0]);
-  output[10] = 2 * (((limb) ((s32) in2[5])) * ((s32) in[5]) +
-                    ((limb) ((s32) in2[3])) * ((s32) in[7]) +
-                    ((limb) ((s32) in2[7])) * ((s32) in[3]) +
-                    ((limb) ((s32) in2[1])) * ((s32) in[9]) +
-                    ((limb) ((s32) in2[9])) * ((s32) in[1])) +
-                    ((limb) ((s32) in2[4])) * ((s32) in[6]) +
-                    ((limb) ((s32) in2[6])) * ((s32) in[4]) +
-                    ((limb) ((s32) in2[2])) * ((s32) in[8]) +
-                    ((limb) ((s32) in2[8])) * ((s32) in[2]);
-  output[11] =      ((limb) ((s32) in2[5])) * ((s32) in[6]) +
-                    ((limb) ((s32) in2[6])) * ((s32) in[5]) +
-                    ((limb) ((s32) in2[4])) * ((s32) in[7]) +
-                    ((limb) ((s32) in2[7])) * ((s32) in[4]) +
-                    ((limb) ((s32) in2[3])) * ((s32) in[8]) +
-                    ((limb) ((s32) in2[8])) * ((s32) in[3]) +
-                    ((limb) ((s32) in2[2])) * ((s32) in[9]) +
-                    ((limb) ((s32) in2[9])) * ((s32) in[2]);
-  output[12] =      ((limb) ((s32) in2[6])) * ((s32) in[6]) +
-               2 * (((limb) ((s32) in2[5])) * ((s32) in[7]) +
-                    ((limb) ((s32) in2[7])) * ((s32) in[5]) +
-                    ((limb) ((s32) in2[3])) * ((s32) in[9]) +
-                    ((limb) ((s32) in2[9])) * ((s32) in[3])) +
-                    ((limb) ((s32) in2[4])) * ((s32) in[8]) +
-                    ((limb) ((s32) in2[8])) * ((s32) in[4]);
-  output[13] =      ((limb) ((s32) in2[6])) * ((s32) in[7]) +
-                    ((limb) ((s32) in2[7])) * ((s32) in[6]) +
-                    ((limb) ((s32) in2[5])) * ((s32) in[8]) +
-                    ((limb) ((s32) in2[8])) * ((s32) in[5]) +
-                    ((limb) ((s32) in2[4])) * ((s32) in[9]) +
-                    ((limb) ((s32) in2[9])) * ((s32) in[4]);
-  output[14] = 2 * (((limb) ((s32) in2[7])) * ((s32) in[7]) +
-                    ((limb) ((s32) in2[5])) * ((s32) in[9]) +
-                    ((limb) ((s32) in2[9])) * ((s32) in[5])) +
-                    ((limb) ((s32) in2[6])) * ((s32) in[8]) +
-                    ((limb) ((s32) in2[8])) * ((s32) in[6]);
-  output[15] =      ((limb) ((s32) in2[7])) * ((s32) in[8]) +
-                    ((limb) ((s32) in2[8])) * ((s32) in[7]) +
-                    ((limb) ((s32) in2[6])) * ((s32) in[9]) +
-                    ((limb) ((s32) in2[9])) * ((s32) in[6]);
-  output[16] =      ((limb) ((s32) in2[8])) * ((s32) in[8]) +
-               2 * (((limb) ((s32) in2[7])) * ((s32) in[9]) +
-                    ((limb) ((s32) in2[9])) * ((s32) in[7]));
-  output[17] =      ((limb) ((s32) in2[8])) * ((s32) in[9]) +
-                    ((limb) ((s32) in2[9])) * ((s32) in[8]);
-  output[18] = 2 *  ((limb) ((s32) in2[9])) * ((s32) in[9]);
-}
-
-/* Reduce a long form to a short form by taking the input mod 2^255 - 19. */
-static void freduce_degree(limb *output) {
-  /* Each of these shifts and adds ends up multiplying the value by 19. */
-  output[8] += output[18] << 4;
-  output[8] += output[18] << 1;
-  output[8] += output[18];
-  output[7] += output[17] << 4;
-  output[7] += output[17] << 1;
-  output[7] += output[17];
-  output[6] += output[16] << 4;
-  output[6] += output[16] << 1;
-  output[6] += output[16];
-  output[5] += output[15] << 4;
-  output[5] += output[15] << 1;
-  output[5] += output[15];
-  output[4] += output[14] << 4;
-  output[4] += output[14] << 1;
-  output[4] += output[14];
-  output[3] += output[13] << 4;
-  output[3] += output[13] << 1;
-  output[3] += output[13];
-  output[2] += output[12] << 4;
-  output[2] += output[12] << 1;
-  output[2] += output[12];
-  output[1] += output[11] << 4;
-  output[1] += output[11] << 1;
-  output[1] += output[11];
-  output[0] += output[10] << 4;
-  output[0] += output[10] << 1;
-  output[0] += output[10];
-}
-
-#if (-1 & 3) != 3
-#error "This code only works on a two's complement system"
-#endif
-
-/* return v / 2^26, using only shifts and adds. */
-static limb div_by_2_26(const limb v)
-{
-  /* High word of v; no shift needed*/
-  const uint32_t highword = (uint32_t) (((uint64_t) v) >> 32);
-  /* Set to all 1s if v was negative; else set to 0s. */
-  const int32_t sign = ((int32_t) highword) >> 31;
-  /* Set to 0x3ffffff if v was negative; else set to 0. */
-  const int32_t roundoff = ((uint32_t) sign) >> 6;
-  /* Should return v / (1<<26) */
-  return (v + roundoff) >> 26;
-}
-
-/* return v / (2^25), using only shifts and adds. */
-static limb div_by_2_25(const limb v)
-{
-  /* High word of v; no shift needed*/
-  const uint32_t highword = (uint32_t) (((uint64_t) v) >> 32);
-  /* Set to all 1s if v was negative; else set to 0s. */
-  const int32_t sign = ((int32_t) highword) >> 31;
-  /* Set to 0x1ffffff if v was negative; else set to 0. */
-  const int32_t roundoff = ((uint32_t) sign) >> 7;
-  /* Should return v / (1<<25) */
-  return (v + roundoff) >> 25;
-}
-
-static s32 div_s32_by_2_25(const s32 v)
-{
-   const s32 roundoff = ((uint32_t)(v >> 31)) >> 7;
-   return (v + roundoff) >> 25;
-}
-
-/* Reduce all coefficients of the short form input so that |x| < 2^26.
- *
- * On entry: |output[i]| < 2^62
- */
-static void freduce_coefficients(limb *output) {
-  unsigned i;
-
-  output[10] = 0;
-
-  for (i = 0; i < 10; i += 2) {
-    limb over = div_by_2_26(output[i]);
-    output[i] -= over << 26;
-    output[i+1] += over;
-
-    over = div_by_2_25(output[i+1]);
-    output[i+1] -= over << 25;
-    output[i+2] += over;
-  }
-  /* Now |output[10]| < 2 ^ 38 and all other coefficients are reduced. */
-  output[0] += output[10] << 4;
-  output[0] += output[10] << 1;
-  output[0] += output[10];
-
-  output[10] = 0;
-
-  /* Now output[1..9] are reduced, and |output[0]| < 2^26 + 19 * 2^38
-   * So |over| will be no more than 77825  */
-  {
-    limb over = div_by_2_26(output[0]);
-    output[0] -= over << 26;
-    output[1] += over;
-  }
-
-  /* Now output[0,2..9] are reduced, and |output[1]| < 2^25 + 77825
-   * So |over| will be no more than 1. */
-  {
-    /* output[1] fits in 32 bits, so we can use div_s32_by_2_25 here. */
-    s32 over32 = div_s32_by_2_25((s32) output[1]);
-    output[1] -= over32 << 25;
-    output[2] += over32;
-  }
-
-  /* Finally, output[0,1,3..9] are reduced, and output[2] is "nearly reduced":
-   * we have |output[2]| <= 2^26.  This is good enough for all of our math,
-   * but it will require an extra freduce_coefficients before fcontract. */
-}
-
-/* A helpful wrapper around fproduct: output = in * in2.
- *
- * output must be distinct to both inputs. The output is reduced degree and
- * reduced coefficient.
- */
-static void
-fmul(limb *output, const limb *in, const limb *in2) {
-  limb t[19];
-  fproduct(t, in, in2);
-  freduce_degree(t);
-  freduce_coefficients(t);
-  memcpy(output, t, sizeof(limb) * 10);
-}
-
-static void fsquare_inner(limb *output, const limb *in) {
-  output[0] =       ((limb) ((s32) in[0])) * ((s32) in[0]);
-  output[1] =  2 *  ((limb) ((s32) in[0])) * ((s32) in[1]);
-  output[2] =  2 * (((limb) ((s32) in[1])) * ((s32) in[1]) +
-                    ((limb) ((s32) in[0])) * ((s32) in[2]));
-  output[3] =  2 * (((limb) ((s32) in[1])) * ((s32) in[2]) +
-                    ((limb) ((s32) in[0])) * ((s32) in[3]));
-  output[4] =       ((limb) ((s32) in[2])) * ((s32) in[2]) +
-               4 *  ((limb) ((s32) in[1])) * ((s32) in[3]) +
-               2 *  ((limb) ((s32) in[0])) * ((s32) in[4]);
-  output[5] =  2 * (((limb) ((s32) in[2])) * ((s32) in[3]) +
-                    ((limb) ((s32) in[1])) * ((s32) in[4]) +
-                    ((limb) ((s32) in[0])) * ((s32) in[5]));
-  output[6] =  2 * (((limb) ((s32) in[3])) * ((s32) in[3]) +
-                    ((limb) ((s32) in[2])) * ((s32) in[4]) +
-                    ((limb) ((s32) in[0])) * ((s32) in[6]) +
-               2 *  ((limb) ((s32) in[1])) * ((s32) in[5]));
-  output[7] =  2 * (((limb) ((s32) in[3])) * ((s32) in[4]) +
-                    ((limb) ((s32) in[2])) * ((s32) in[5]) +
-                    ((limb) ((s32) in[1])) * ((s32) in[6]) +
-                    ((limb) ((s32) in[0])) * ((s32) in[7]));
-  output[8] =       ((limb) ((s32) in[4])) * ((s32) in[4]) +
-               2 * (((limb) ((s32) in[2])) * ((s32) in[6]) +
-                    ((limb) ((s32) in[0])) * ((s32) in[8]) +
-               2 * (((limb) ((s32) in[1])) * ((s32) in[7]) +
-                    ((limb) ((s32) in[3])) * ((s32) in[5])));
-  output[9] =  2 * (((limb) ((s32) in[4])) * ((s32) in[5]) +
-                    ((limb) ((s32) in[3])) * ((s32) in[6]) +
-                    ((limb) ((s32) in[2])) * ((s32) in[7]) +
-                    ((limb) ((s32) in[1])) * ((s32) in[8]) +
-                    ((limb) ((s32) in[0])) * ((s32) in[9]));
-  output[10] = 2 * (((limb) ((s32) in[5])) * ((s32) in[5]) +
-                    ((limb) ((s32) in[4])) * ((s32) in[6]) +
-                    ((limb) ((s32) in[2])) * ((s32) in[8]) +
-               2 * (((limb) ((s32) in[3])) * ((s32) in[7]) +
-                    ((limb) ((s32) in[1])) * ((s32) in[9])));
-  output[11] = 2 * (((limb) ((s32) in[5])) * ((s32) in[6]) +
-                    ((limb) ((s32) in[4])) * ((s32) in[7]) +
-                    ((limb) ((s32) in[3])) * ((s32) in[8]) +
-                    ((limb) ((s32) in[2])) * ((s32) in[9]));
-  output[12] =      ((limb) ((s32) in[6])) * ((s32) in[6]) +
-               2 * (((limb) ((s32) in[4])) * ((s32) in[8]) +
-               2 * (((limb) ((s32) in[5])) * ((s32) in[7]) +
-                    ((limb) ((s32) in[3])) * ((s32) in[9])));
-  output[13] = 2 * (((limb) ((s32) in[6])) * ((s32) in[7]) +
-                    ((limb) ((s32) in[5])) * ((s32) in[8]) +
-                    ((limb) ((s32) in[4])) * ((s32) in[9]));
-  output[14] = 2 * (((limb) ((s32) in[7])) * ((s32) in[7]) +
-                    ((limb) ((s32) in[6])) * ((s32) in[8]) +
-               2 *  ((limb) ((s32) in[5])) * ((s32) in[9]));
-  output[15] = 2 * (((limb) ((s32) in[7])) * ((s32) in[8]) +
-                    ((limb) ((s32) in[6])) * ((s32) in[9]));
-  output[16] =      ((limb) ((s32) in[8])) * ((s32) in[8]) +
-               4 *  ((limb) ((s32) in[7])) * ((s32) in[9]);
-  output[17] = 2 *  ((limb) ((s32) in[8])) * ((s32) in[9]);
-  output[18] = 2 *  ((limb) ((s32) in[9])) * ((s32) in[9]);
-}
-
-static void
-fsquare(limb *output, const limb *in) {
-  limb t[19];
-  fsquare_inner(t, in);
-  freduce_degree(t);
-  freduce_coefficients(t);
-  memcpy(output, t, sizeof(limb) * 10);
-}
-
-/* Take a little-endian, 32-byte number and expand it into polynomial form */
-static void
-fexpand(limb *output, const u8 *input) {
-#define F(n,start,shift,mask) \
-  output[n] = ((((limb) input[start + 0]) | \
-                ((limb) input[start + 1]) << 8 | \
-                ((limb) input[start + 2]) << 16 | \
-                ((limb) input[start + 3]) << 24) >> shift) & mask;
-  F(0, 0, 0, 0x3ffffff);
-  F(1, 3, 2, 0x1ffffff);
-  F(2, 6, 3, 0x3ffffff);
-  F(3, 9, 5, 0x1ffffff);
-  F(4, 12, 6, 0x3ffffff);
-  F(5, 16, 0, 0x1ffffff);
-  F(6, 19, 1, 0x3ffffff);
-  F(7, 22, 3, 0x1ffffff);
-  F(8, 25, 4, 0x3ffffff);
-  F(9, 28, 6, 0x1ffffff);
-#undef F
-}
-
-#if (-32 >> 1) != -16
-#error "This code only works when >> does sign-extension on negative numbers"
-#endif
-
-/* Take a fully reduced polynomial form number and contract it into a
- * little-endian, 32-byte array
- */
-static void
-fcontract(u8 *output, limb *input) {
-  int i;
-  int j;
-
-  for (j = 0; j < 2; ++j) {
-    for (i = 0; i < 9; ++i) {
-      if ((i & 1) == 1) {
-        /* This calculation is a time-invariant way to make input[i] positive
-           by borrowing from the next-larger limb.
-        */
-        const s32 mask = (s32)(input[i]) >> 31;
-        const s32 carry = -(((s32)(input[i]) & mask) >> 25);
-        input[i] = (s32)(input[i]) + (carry << 25);
-        input[i+1] = (s32)(input[i+1]) - carry;
-      } else {
-        const s32 mask = (s32)(input[i]) >> 31;
-        const s32 carry = -(((s32)(input[i]) & mask) >> 26);
-        input[i] = (s32)(input[i]) + (carry << 26);
-        input[i+1] = (s32)(input[i+1]) - carry;
-      }
-    }
-    {
-      const s32 mask = (s32)(input[9]) >> 31;
-      const s32 carry = -(((s32)(input[9]) & mask) >> 25);
-      input[9] = (s32)(input[9]) + (carry << 25);
-      input[0] = (s32)(input[0]) - (carry * 19);
-    }
-  }
-
-  /* The first borrow-propagation pass above ended with every limb
-     except (possibly) input[0] non-negative.
-
-     Since each input limb except input[0] is decreased by at most 1
-     by a borrow-propagation pass, the second borrow-propagation pass
-     could only have wrapped around to decrease input[0] again if the
-     first pass left input[0] negative *and* input[1] through input[9]
-     were all zero.  In that case, input[1] is now 2^25 - 1, and this
-     last borrow-propagation step will leave input[1] non-negative.
-  */
-  {
-    const s32 mask = (s32)(input[0]) >> 31;
-    const s32 carry = -(((s32)(input[0]) & mask) >> 26);
-    input[0] = (s32)(input[0]) + (carry << 26);
-    input[1] = (s32)(input[1]) - carry;
-  }
-
-  /* Both passes through the above loop, plus the last 0-to-1 step, are
-     necessary: if input[9] is -1 and input[0] through input[8] are 0,
-     negative values will remain in the array until the end.
-   */
-
-  input[1] <<= 2;
-  input[2] <<= 3;
-  input[3] <<= 5;
-  input[4] <<= 6;
-  input[6] <<= 1;
-  input[7] <<= 3;
-  input[8] <<= 4;
-  input[9] <<= 6;
-#define F(i, s) \
-  output[s+0] |=  input[i] & 0xff; \
-  output[s+1]  = (input[i] >> 8) & 0xff; \
-  output[s+2]  = (input[i] >> 16) & 0xff; \
-  output[s+3]  = (input[i] >> 24) & 0xff;
-  output[0] = 0;
-  output[16] = 0;
-  F(0,0);
-  F(1,3);
-  F(2,6);
-  F(3,9);
-  F(4,12);
-  F(5,16);
-  F(6,19);
-  F(7,22);
-  F(8,25);
-  F(9,28);
-#undef F
-}
-
-/* Input: Q, Q', Q-Q'
- * Output: 2Q, Q+Q'
- *
- *   x2 z3: long form
- *   x3 z3: long form
- *   x z: short form, destroyed
- *   xprime zprime: short form, destroyed
- *   qmqp: short form, preserved
- */
-static void fmonty(limb *x2, limb *z2,  /* output 2Q */
-                   limb *x3, limb *z3,  /* output Q + Q' */
-                   limb *x, limb *z,    /* input Q */
-                   limb *xprime, limb *zprime,  /* input Q' */
-                   const limb *qmqp /* input Q - Q' */) {
-  limb origx[10], origxprime[10], zzz[19], xx[19], zz[19], xxprime[19],
-        zzprime[19], zzzprime[19], xxxprime[19];
-
-  memcpy(origx, x, 10 * sizeof(limb));
-  fsum(x, z);
-  fdifference(z, origx);  /* does x - z */
-
-  memcpy(origxprime, xprime, sizeof(limb) * 10);
-  fsum(xprime, zprime);
-  fdifference(zprime, origxprime);
-  fproduct(xxprime, xprime, z);
-  fproduct(zzprime, x, zprime);
-  freduce_degree(xxprime);
-  freduce_coefficients(xxprime);
-  freduce_degree(zzprime);
-  freduce_coefficients(zzprime);
-  memcpy(origxprime, xxprime, sizeof(limb) * 10);
-  fsum(xxprime, zzprime);
-  fdifference(zzprime, origxprime);
-  fsquare(xxxprime, xxprime);
-  fsquare(zzzprime, zzprime);
-  fproduct(zzprime, zzzprime, qmqp);
-  freduce_degree(zzprime);
-  freduce_coefficients(zzprime);
-  memcpy(x3, xxxprime, sizeof(limb) * 10);
-  memcpy(z3, zzprime, sizeof(limb) * 10);
-
-  fsquare(xx, x);
-  fsquare(zz, z);
-  fproduct(x2, xx, zz);
-  freduce_degree(x2);
-  freduce_coefficients(x2);
-  fdifference(zz, xx);         /* does zz = xx - zz */
-  memset(zzz + 10, 0, sizeof(limb) * 9);
-  fscalar_product(zzz, zz, 121665);
-  /* No need to call freduce_degree here:
-     fscalar_product doesn't increase the degree of its input.
-   */
-  freduce_coefficients(zzz);
-  fsum(zzz, xx);
-  fproduct(z2, zz, zzz);
-  freduce_degree(z2);
-  freduce_coefficients(z2);
-}
-
-/* Conditionally swap two reduced-form limb arrays if 'iswap' is 1, but leave
- * them unchanged if 'iswap' is 0.  Runs in data-invariant time to avoid
- * side-channel attacks.
- *
- * NOTE that this function requires that 'iswap' be 1 or 0; other values give
- * wrong results.  Also, the two limb arrays must be in reduced-coefficient,
- * reduced-degree form: the values in a[10..19] or b[10..19] aren't swapped,
- * and all all values in a[0..9],b[0..9] must have magnitude less than
- * INT32_MAX.
- */
-static void
-swap_conditional(limb a[19], limb b[19], limb iswap) {
-  unsigned i;
-  const s32 swap = (s32) -iswap;
-
-  for (i = 0; i < 10; ++i) {
-    const s32 x = swap & ( ((s32)a[i]) ^ ((s32)b[i]) );
-    a[i] = ((s32)a[i]) ^ x;
-    b[i] = ((s32)b[i]) ^ x;
-  }
-}
-
-/* Calculates nQ where Q is the x-coordinate of a point on the curve
- *
- *   resultx/resultz: the x coordinate of the resulting curve point (short form)
- *   n: a little endian, 32-byte number
- *   q: a point of the curve (short form)
- */
-static void
-cmult(limb *resultx, limb *resultz, const u8 *n, const limb *q) {
-  limb a[19] = {0}, b[19] = {1}, c[19] = {1}, d[19] = {0};
-  limb *nqpqx = a, *nqpqz = b, *nqx = c, *nqz = d, *t;
-  limb e[19] = {0}, f[19] = {1}, g[19] = {0}, h[19] = {1};
-  limb *nqpqx2 = e, *nqpqz2 = f, *nqx2 = g, *nqz2 = h;
-
-  unsigned i, j;
-
-  memcpy(nqpqx, q, sizeof(limb) * 10);
-
-  for (i = 0; i < 32; ++i) {
-    u8 byte = n[31 - i];
-    for (j = 0; j < 8; ++j) {
-      const limb bit = byte >> 7;
-
-      swap_conditional(nqx, nqpqx, bit);
-      swap_conditional(nqz, nqpqz, bit);
-      fmonty(nqx2, nqz2,
-             nqpqx2, nqpqz2,
-             nqx, nqz,
-             nqpqx, nqpqz,
-             q);
-      swap_conditional(nqx2, nqpqx2, bit);
-      swap_conditional(nqz2, nqpqz2, bit);
-
-      t = nqx;
-      nqx = nqx2;
-      nqx2 = t;
-      t = nqz;
-      nqz = nqz2;
-      nqz2 = t;
-      t = nqpqx;
-      nqpqx = nqpqx2;
-      nqpqx2 = t;
-      t = nqpqz;
-      nqpqz = nqpqz2;
-      nqpqz2 = t;
-
-      byte <<= 1;
-    }
-  }
-
-  memcpy(resultx, nqx, sizeof(limb) * 10);
-  memcpy(resultz, nqz, sizeof(limb) * 10);
-}
-
-/* -----------------------------------------------------------------------------
- * Shamelessly copied from djb's code
- * ----------------------------------------------------------------------------- */
-static void
-crecip(limb *out, const limb *z) {
-  limb z2[10];
-  limb z9[10];
-  limb z11[10];
-  limb z2_5_0[10];
-  limb z2_10_0[10];
-  limb z2_20_0[10];
-  limb z2_50_0[10];
-  limb z2_100_0[10];
-  limb t0[10];
-  limb t1[10];
-  int i;
-
-  /* 2 */ fsquare(z2,z);
-  /* 4 */ fsquare(t1,z2);
-  /* 8 */ fsquare(t0,t1);
-  /* 9 */ fmul(z9,t0,z);
-  /* 11 */ fmul(z11,z9,z2);
-  /* 22 */ fsquare(t0,z11);
-  /* 2^5 - 2^0 = 31 */ fmul(z2_5_0,t0,z9);
-
-  /* 2^6 - 2^1 */ fsquare(t0,z2_5_0);
-  /* 2^7 - 2^2 */ fsquare(t1,t0);
-  /* 2^8 - 2^3 */ fsquare(t0,t1);
-  /* 2^9 - 2^4 */ fsquare(t1,t0);
-  /* 2^10 - 2^5 */ fsquare(t0,t1);
-  /* 2^10 - 2^0 */ fmul(z2_10_0,t0,z2_5_0);
-
-  /* 2^11 - 2^1 */ fsquare(t0,z2_10_0);
-  /* 2^12 - 2^2 */ fsquare(t1,t0);
-  /* 2^20 - 2^10 */ for (i = 2;i < 10;i += 2) { fsquare(t0,t1); fsquare(t1,t0); }
-  /* 2^20 - 2^0 */ fmul(z2_20_0,t1,z2_10_0);
-
-  /* 2^21 - 2^1 */ fsquare(t0,z2_20_0);
-  /* 2^22 - 2^2 */ fsquare(t1,t0);
-  /* 2^40 - 2^20 */ for (i = 2;i < 20;i += 2) { fsquare(t0,t1); fsquare(t1,t0); }
-  /* 2^40 - 2^0 */ fmul(t0,t1,z2_20_0);
-
-  /* 2^41 - 2^1 */ fsquare(t1,t0);
-  /* 2^42 - 2^2 */ fsquare(t0,t1);
-  /* 2^50 - 2^10 */ for (i = 2;i < 10;i += 2) { fsquare(t1,t0); fsquare(t0,t1); }
-  /* 2^50 - 2^0 */ fmul(z2_50_0,t0,z2_10_0);
-
-  /* 2^51 - 2^1 */ fsquare(t0,z2_50_0);
-  /* 2^52 - 2^2 */ fsquare(t1,t0);
-  /* 2^100 - 2^50 */ for (i = 2;i < 50;i += 2) { fsquare(t0,t1); fsquare(t1,t0); }
-  /* 2^100 - 2^0 */ fmul(z2_100_0,t1,z2_50_0);
-
-  /* 2^101 - 2^1 */ fsquare(t1,z2_100_0);
-  /* 2^102 - 2^2 */ fsquare(t0,t1);
-  /* 2^200 - 2^100 */ for (i = 2;i < 100;i += 2) { fsquare(t1,t0); fsquare(t0,t1); }
-  /* 2^200 - 2^0 */ fmul(t1,t0,z2_100_0);
-
-  /* 2^201 - 2^1 */ fsquare(t0,t1);
-  /* 2^202 - 2^2 */ fsquare(t1,t0);
-  /* 2^250 - 2^50 */ for (i = 2;i < 50;i += 2) { fsquare(t0,t1); fsquare(t1,t0); }
-  /* 2^250 - 2^0 */ fmul(t0,t1,z2_50_0);
-
-  /* 2^251 - 2^1 */ fsquare(t1,t0);
-  /* 2^252 - 2^2 */ fsquare(t0,t1);
-  /* 2^253 - 2^3 */ fsquare(t1,t0);
-  /* 2^254 - 2^4 */ fsquare(t0,t1);
-  /* 2^255 - 2^5 */ fsquare(t1,t0);
-  /* 2^255 - 21 */ fmul(out,t1,z11);
-}
-
-int curve25519_donna(u8 *, const u8 *, const u8 *);
-
-int curve25519_donna(u8 *mypublic, const u8 *secret, const u8 *basepoint) {
-  limb bp[10], x[10], z[11], zmone[10];
-  uint8_t e[32];
-  int i;
-
-  for (i = 0; i < 32; ++i) e[i] = secret[i];
-  e[0] &= 248;
-  e[31] &= 127;
-  e[31] |= 64;
-
-  fexpand(bp, basepoint);
-  cmult(x, z, e, bp);
-  crecip(zmone, z);
-  fmul(z, x, zmone);
-  freduce_coefficients(z);
-  fcontract(mypublic, z);
-  return 0;
-}
diff --git a/jni/libzrtp/sources/bnlib/ec/ec.c b/jni/libzrtp/sources/bnlib/ec/ec.c
index 18e612f..1f4123a 100644
--- a/jni/libzrtp/sources/bnlib/ec/ec.c
+++ b/jni/libzrtp/sources/bnlib/ec/ec.c
@@ -1,5 +1,5 @@
 /*
- * Copyright (C) 2012-2013 Werner Dittmann
+ * Copyright (C) 2012 Werner Dittmann
  * All rights reserved. For licensing and other legal details, see the file legal.c.
  *
  * @author Werner Dittmann <Werner.Dittmann@t-online.de>
@@ -100,46 +100,6 @@
         "11839296a789a3bc0045c8a5fb42c7d1bd998f54449579b446817afbd17273e662c97ee72995ef42640c550b9013fad0761353c7086a272c24088be94769fd16650",
 };
 
-
-/*
- * The data for curve3617 copied from:
- * http://safecurves.cr.yp.to/field.html
- * http://safecurves.cr.yp.to/base.html
- */
-static curveData curve3617 = {
-    "3fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffef",  /* Prime */
-    "7ffffffffffffffffffffffffffffffffffffffffffffffffffeb3cc92414cf706022b36f1c0338ad63cf181b0e71a5e106af79",   /* order */
-    "",                                                                                                          /* SEED */
-    "",                                                                                                          /* c */
-    "",                                                                                                          /* b */
-    "1a334905141443300218c0631c326e5fcd46369f44c03ec7f57ff35498a4ab4d6d6ba111301a73faa8537c64c4fd3812f3cbc595",  /* Gx*/
-    "22",                                                                                                        /* Gy (radix 16) */
-};
-
-/*
- * The data for curve25519 copied from:
- * http://safecurves.cr.yp.to/field.html
- * http://safecurves.cr.yp.to/base.html
- * 
- * Note: 
- * The data for Curve25519 is here for the sake of completeness and to have the same
- * set of initialization. One exception if the base point X coordinate (Gx) that we use to
- * compute the DH public value, refer to function ecdhGeneratePublic(...) in ecdh.c.
- * 
- * Otherwise the functions use EcCurve structure only to get the pointers to the Curve25519
- * wrapper functions.
- * 
- */
-static curveData curve25519 = {
-    "7fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffed",   /* Prime */
-    "1000000000000000000000000000000014def9dea2f79cd65812631a5cf5d3ed",   /* order */
-    "",                                                                   /* SEED */
-    "",                                                                   /* c */
-    "",                                                                   /* b */
-    "9",                                                                  /* Gx */
-    "20ae19a1b8a086b4e01edd2c7748d14c923d4d7e6d7c61b229e9c5a27eced3d9",   /* Gy */
-};
-
 /*============================================================================*/
 /*    Bignum Shorthand Functions                                              */
 /*============================================================================*/
@@ -180,90 +140,53 @@
     return 0;
 }
 
-int bnMulMod_ (struct BigNum *rslt, struct BigNum *n1, struct BigNum *n2, struct BigNum *mod, const EcCurve *curve)
+int bnMulMod_ (struct BigNum *rslt, struct BigNum *n1, struct BigNum *n2, struct BigNum *mod)
 {
     bnMul (rslt, n1, n2);
-    if (curve)
-        curve->modOp(rslt, rslt, mod);
-    else
-        bnMod(rslt, rslt, mod);
+    bnMod (rslt, rslt, mod);
     return 0;
 }
 
-int bnMulQMod_ (struct BigNum *rslt, struct BigNum *n1, unsigned n2, struct BigNum *mod, const EcCurve *curve)
+int bnMulQMod_ (struct BigNum *rslt, struct BigNum *n1, unsigned n2, struct BigNum *mod)
 {
     bnMulQ (rslt, n1, n2);
-    if (curve)
-        curve->modOp(rslt, rslt, mod);
-    else
-        bnMod(rslt, rslt, mod);
-   return 0;
-}
-
-int bnSquareMod_ (struct BigNum *rslt, struct BigNum *n1, struct BigNum *mod, const EcCurve *curve)
-{
-    bnSquare (rslt, n1);
-    if (curve)
-        curve->modOp(rslt, rslt, mod);
-    else
-        bnMod(rslt, rslt, mod);
+    bnMod (rslt, rslt, mod);
     return 0;
 }
 
-/*
- * Note on the Curve25519 functions and usage of BigNumber:
- * In most cases the functions to compute Curve25519 data are small wrapper functions
- * that implement the same API as for the other curve functions. The wrapper functions
- * then call the very specific, high-efficient function in curve25519-donna.c .
- * 
- * For Curve25519 we don't have a real implementation for point add, point doubling, modulo
- * and check public key. Please refer to the actual implementations below.
- */
-
-static int ecGetAffineNist(const EcCurve *curve, EcPoint *R, const EcPoint *P);
-static int ecGetAffineEd(const EcCurve *curve, EcPoint *R, const EcPoint *P);
-static int ecGetAffine25519(const EcCurve *curve, EcPoint *R, const EcPoint *P);
-
-static int ecDoublePointNist(const EcCurve *curve, EcPoint *R, const EcPoint *P);
-static int ecDoublePointEd(const EcCurve *curve, EcPoint *R, const EcPoint *P);
-static int ecDoublePoint25519(const EcCurve *curve, EcPoint *R, const EcPoint *P);
-
-static int ecAddPointNist(const EcCurve *curve, EcPoint *R, const EcPoint *P, const EcPoint *Q);
-static int ecAddPointEd(const EcCurve *curve, EcPoint *R, const EcPoint *P, const EcPoint *Q);
-static int ecAddPoint25519(const EcCurve *curve, EcPoint *R, const EcPoint *P, const EcPoint *Q);
-
-static int ecCheckPubKeyNist(const EcCurve *curve, const EcPoint *pub);
-static int ecCheckPubKey3617(const EcCurve *curve, const EcPoint *pub);
-static int ecCheckPubKey25519(const EcCurve *curve, const EcPoint *pub);
-
-static int ecGenerateRandomNumberNist(const EcCurve *curve, BigNum *d);
-static int ecGenerateRandomNumber3617(const EcCurve *curve, BigNum *d);
-static int ecGenerateRandomNumber25519(const EcCurve *curve, BigNum *d);
-
-static int ecMulPointScalarNormal(const EcCurve *curve, EcPoint *R, const EcPoint *P, const BigNum *scalar);
-static int ecMulPointScalar25519(const EcCurve *curve, EcPoint *R, const EcPoint *P, const BigNum *scalar);
-
-/* Forward declaration of new modulo functions for the EC curves */
-static int newMod192(BigNum *r, const BigNum *a, const BigNum *modulo);
-static int newMod256(BigNum *r, const BigNum *a, const BigNum *modulo);
-static int newMod384(BigNum *r, const BigNum *a, const BigNum *modulo);
-static int newMod521(BigNum *r, const BigNum *a, const BigNum *modulo);
-
-static int mod3617(BigNum *r, const BigNum *a, const BigNum *modulo);
-static int mod25519(BigNum *r, const BigNum *a, const BigNum *modulo);
-
-static void commonInit()
+int bnSquareMod_ (struct BigNum *rslt, struct BigNum *n1, struct BigNum *mod)
 {
-    bnBegin(mpiZero); bnSetQ(mpiZero, 0);
-    bnBegin(mpiOne); bnSetQ(mpiOne, 1);
-    bnBegin(mpiTwo); bnSetQ(mpiTwo, 2);
-    bnBegin(mpiThree); bnSetQ(mpiThree, 3);
-    bnBegin(mpiFour); bnSetQ(mpiFour, 4);
-    bnBegin(mpiEight); bnSetQ(mpiEight, 8);
+    bnSquare (rslt, n1);
+    bnMod (rslt, rslt, mod);
+    return 0;
 }
 
-static void curveCommonInit(EcCurve *curve)
+int ecGetCurveNistECp(NistCurves curveId, NistECpCurve *curve)
 {
+    size_t maxBits;
+    curveData *cd;
+
+    if (!initialized) {
+        bnBegin(mpiZero); bnSetQ(mpiZero, 0);
+        bnBegin(mpiOne); bnSetQ(mpiOne, 1);
+        bnBegin(mpiTwo); bnSetQ(mpiTwo, 2);
+        bnBegin(mpiThree); bnSetQ(mpiThree, 3);
+        bnBegin(mpiFour); bnSetQ(mpiFour, 4);
+        bnBegin(mpiEight); bnSetQ(mpiEight, 8);
+        initialized = 1;
+    }
+    if (curve == NULL)
+        return -2;
+
+    bnBegin(&curve->_p);    curve->p = &curve->_p;
+    bnBegin(&curve->_n);    curve->n = &curve->_n;
+    bnBegin(&curve->_SEED); curve->SEED = &curve->_SEED;
+    bnBegin(&curve->_c);    curve->c = &curve->_c;
+    bnBegin(&curve->_a);    curve->a = &curve->_a;
+    bnBegin(&curve->_b);    curve->b = &curve->_b;
+    bnBegin(&curve->_Gx);   curve->Gx = &curve->_Gx;
+    bnBegin(&curve->_Gy);   curve->Gy = &curve->_Gy;
+
     /* Initialize scratchpad variables and their pointers */
     bnBegin(&curve->_S1); curve->S1 = &curve->_S1;
     bnBegin(&curve->_U1); curve->U1 = &curve->_U1;
@@ -273,11 +196,41 @@
     bnBegin(&curve->_t1); curve->t1 = &curve->_t1;
     bnBegin(&curve->_t2); curve->t2 = &curve->_t2;
     bnBegin(&curve->_t3); curve->t3 = &curve->_t3;
-}
 
-static void curveCommonPrealloc(EcCurve *curve)
-{
-    size_t maxBits;
+    switch (curveId) {
+    case NIST192P:
+        cd = &nist192;
+        break;
+
+    case NIST224P:
+        cd = &nist224;
+        break;
+
+    case NIST256P:
+        cd = &nist256;
+        break;
+
+    case NIST384P:
+        cd = &nist384;
+        break;
+
+    case NIST521P:
+        cd = &nist521;
+        break;
+
+    default:
+        return -2;
+    }
+
+    bnReadAscii(curve->p, cd->p, 10);
+    bnReadAscii(curve->n, cd->n, 10);
+    bnReadAscii(curve->SEED, cd->SEED, 16);
+    bnReadAscii(curve->c, cd->c, 16);
+    bnCopy(curve->a, curve->p);
+    bnSub(curve->a, mpiThree);
+    bnReadAscii(curve->b, cd->b, 16);
+    bnReadAscii(curve->Gx, cd->Gx, 16);
+    bnReadAscii(curve->Gy, cd->Gy, 16);
 
     /* variables must be able to hold p^2, plus one nimb (min. 15 bits) for overflow */
     maxBits = bnBits(curve->p) * 2 + 15;
@@ -292,151 +245,16 @@
     bnPrealloc(curve->t1, maxBits);
     bnPrealloc(curve->t2, maxBits);
     bnPrealloc(curve->t3, maxBits);
-}
-
-int ecGetCurveNistECp(Curves curveId, EcCurve *curve)
-{
-    curveData *cd;
-
-    if (curveId >= Curve25519 && curveId <= Curve3617)
-        return ecGetCurvesCurve(curveId, curve);
-
-    if (!initialized) {
-        commonInit();
-        initialized = 1;
-    }
-    if (curve == NULL)
-        return -2;
-
-    bnBegin(&curve->_p);    curve->p = &curve->_p;
-    bnBegin(&curve->_n);    curve->n = &curve->_n;
-    bnBegin(&curve->_SEED); curve->SEED = &curve->_SEED;
-    bnBegin(&curve->_c);    curve->c = &curve->_c;
-    bnBegin(&curve->_a);    curve->a = &curve->_a;
-    bnBegin(&curve->_b);    curve->b = &curve->_b;
-    bnBegin(&curve->_Gx);   curve->Gx = &curve->_Gx;
-    bnBegin(&curve->_Gy);   curve->Gy = &curve->_Gy;
-
-    curveCommonInit(curve);
-
-    switch (curveId) {
-    case NIST192P:
-        cd = &nist192;
-        curve->modOp = newMod192;
-        break;
-
-    case NIST224P:
-        cd = &nist224;
-        curve->modOp = bnMod;
-        break;
-
-    case NIST256P:
-        cd = &nist256;
-        curve->modOp = bnMod;
-        break;
-
-    case NIST384P:
-        cd = &nist384;
-        curve->modOp = newMod384;
-        break;
-
-    case NIST521P:
-        cd = &nist521;
-        curve->modOp = newMod521;
-        break;
-
-    default:
-        return -2;
-    }
-
-    curve->affineOp = ecGetAffineNist;
-    curve->doubleOp = ecDoublePointNist;
-    curve->addOp = ecAddPointNist;
-    curve->checkPubOp = ecCheckPubKeyNist;
-    curve->randomOp = ecGenerateRandomNumberNist;
-    curve->mulScalar = ecMulPointScalarNormal;
-
-    bnReadAscii(curve->p, cd->p, 10);
-    bnReadAscii(curve->n, cd->n, 10);
-    bnReadAscii(curve->SEED, cd->SEED, 16);
-    bnReadAscii(curve->c, cd->c, 16);
-    bnCopy(curve->a, curve->p);
-    bnSub(curve->a, mpiThree);
-    bnReadAscii(curve->b, cd->b, 16);
-    bnReadAscii(curve->Gx, cd->Gx, 16);
-    bnReadAscii(curve->Gy, cd->Gy, 16);
-
-    curveCommonPrealloc(curve);
-    curve->id = curveId;
 
     return 0;
+
+/*     ecFreeCurveNistECp(curve);
+     return ret;
+*/
 }
 
-int ecGetCurvesCurve(Curves curveId, EcCurve *curve)
-{
-    curveData *cd;
 
-    if (!initialized) {
-        commonInit();
-        initialized = 1;
-    }
-    if (curve == NULL)
-        return -2;
-
-    /* set-up all bignum structures, simplifies "free" handling */
-    bnBegin(&curve->_p);    curve->p = &curve->_p;
-    bnBegin(&curve->_n);    curve->n = &curve->_n;
-    bnBegin(&curve->_SEED); curve->SEED = &curve->_SEED;
-    bnBegin(&curve->_c);    curve->c = &curve->_c;
-    bnBegin(&curve->_a);    curve->a = &curve->_a;
-    bnBegin(&curve->_b);    curve->b = &curve->_b;
-    bnBegin(&curve->_Gx);   curve->Gx = &curve->_Gx;
-    bnBegin(&curve->_Gy);   curve->Gy = &curve->_Gy;
-
-    curveCommonInit(curve);
-
-    switch (curveId) {
-    case Curve3617:
-        cd = &curve3617;
-        curve->modOp = mod3617;
-        curve->affineOp = ecGetAffineEd;
-        curve->doubleOp = ecDoublePointEd;
-        curve->addOp = ecAddPointEd;
-        curve->checkPubOp = ecCheckPubKey3617;
-        curve->randomOp = ecGenerateRandomNumber3617;
-        curve->mulScalar = ecMulPointScalarNormal;
-
-        bnReadAscii(curve->a, "3617", 10);
-        break;
-
-    case Curve25519:
-        cd = &curve25519;
-        curve->modOp = mod25519;
-        curve->affineOp = ecGetAffine25519;
-        curve->doubleOp = ecDoublePoint25519;
-        curve->addOp = ecAddPoint25519;
-        curve->checkPubOp = ecCheckPubKey25519;
-        curve->randomOp = ecGenerateRandomNumber25519;
-        curve->mulScalar = ecMulPointScalar25519;
-
-        bnReadAscii(curve->a, "486662", 10);
-        break;
-
-    default:
-        return -2;
-    }
-    bnReadAscii(curve->p, cd->p, 16);
-    bnReadAscii(curve->n, cd->n, 16);
-
-    bnReadAscii(curve->Gx, cd->Gx, 16);
-    bnReadAscii(curve->Gy, cd->Gy, 16);
-
-    curveCommonPrealloc(curve);
-    curve->id = curveId;
-    return 0;
-}
-
-void ecFreeCurveNistECp(EcCurve *curve) 
+void ecFreeCurveNistECp(NistECpCurve *curve) 
 {
     if (curve == NULL)
         return;
@@ -459,40 +277,12 @@
     bnEnd(curve->t3);
 }
 
-/*
- * EC point helper functions
- */
-
-void ecInitPoint(EcPoint *P)
-{
-    INIT_EC_POINT(P);
-}
-
-void ecFreePoint(EcPoint *P)
-{
-    FREE_EC_POINT(P);
-}
-
-void ecSetBasePoint(EcCurve *C, EcPoint *P)
-{
-    SET_EC_BASE_POINT(C, P);
-}
-
-void ecFreeCurvesCurve(EcCurve *curve)
-{
-    ecFreeCurveNistECp(curve);
-}
 
 /*============================================================================*/
 /*    Elliptic Curve arithmetic                                               */
 /*============================================================================*/
 
-int ecGetAffine(const EcCurve *curve, EcPoint *R, const EcPoint *P)
-{
-    return curve->affineOp(curve, R, P);
-}
-
-static int ecGetAffineNist(const EcCurve *curve, EcPoint *R, const EcPoint *P)
+int ecGetAffine(const NistECpCurve *curve, EcPoint *R, const EcPoint *P)
 {
     int ret = 0;
 
@@ -502,13 +292,13 @@
     bnBegin(&z_2);
 
     /* affine x = X / Z^2 */
-    bnInv (&z_1, P->z, curve->p);                 /* z_1 = Z^(-1) */
-    bnMulMod_(&z_2, &z_1, &z_1, curve->p, curve); /* z_2 = Z^(-2) */
-    bnMulMod_(R->x, P->x, &z_2, curve->p, curve);
-
+    bnInv (&z_1, P->z, curve->p);          /* z_1 = Z^(-1) */
+    bnMulMod_(&z_2, &z_1, &z_1, curve->p); /* z_2 = Z^(-2) */
+    bnMulMod_(R->x, P->x, &z_2, curve->p);
+    
     /* affine y = Y / Z^3 */
-    bnMulMod_(&z_2, &z_2, &z_1, curve->p, curve); /* z_2 = Z^(-3) */
-    bnMulMod_(R->y, P->y, &z_2, curve->p, curve);
+    bnMulMod_(&z_2, &z_2, &z_1, curve->p); /* z_2 = Z^(-3) */
+    bnMulMod_(R->y, P->y, &z_2, curve->p);
 
     bnSetQ(R->z, 1);
 
@@ -517,48 +307,7 @@
     return ret;
 }
 
-static int ecGetAffineEd(const EcCurve *curve, EcPoint *R, const EcPoint *P)
-{
-    int ret = 0;
-
-    struct BigNum z_1;
-
-    bnBegin(&z_1);
-
-    /* affine x = X / Z */
-    bnInv (&z_1, P->z, curve->p);                 /* z_1 = Z^(-1) */
-    bnMulMod_(R->x, P->x, &z_1, curve->p, curve);
-
-    /* affine y = Y / Z */
-    bnMulMod_(R->y, P->y, &z_1, curve->p, curve);
-
-    bnSetQ(R->z, 1);
-
-    bnEnd(&z_1);
-    return ret;
-
-}
-
-/* 
- * If the arguments do not point to the same EcPoint then copy P to result.
- * Curve25519 has no specific GetAffine function, it's all inside curve25519-donna
- */
-static int ecGetAffine25519(const EcCurve *curve, EcPoint *R, const EcPoint *P)
-{
-    if (R != P) {
-        bnCopy(R->x, P->x);
-        bnCopy(R->y, P->y);
-        bnCopy(R->z, P->z);
-    }
-    return 0;
-}
-
-int ecDoublePoint(const EcCurve *curve, EcPoint *R, const EcPoint *P)
-{
-    return curve->doubleOp(curve, R, P);
-}
-
-static int ecDoublePointNist(const EcCurve *curve, EcPoint *R, const EcPoint *P)
+int ecDoublePoint(const NistECpCurve *curve, EcPoint *R, const EcPoint *P)
 {
     int ret = 0;
 
@@ -584,37 +333,37 @@
         ptP = P;
 
     /* S = 4*X*Y^2, save Y^2 in t1 for later use */
-    bnMulMod_(curve->t1, ptP->y, ptP->y, curve->p, curve);       /* t1 = Y^2 */
-    bnMulMod_(curve->t0, ptP->x, mpiFour, curve->p, curve);      /* t0 = 4 * X */
-    bnMulMod_(curve->S1, curve->t0, curve->t1, curve->p, curve); /* S1 = t0 * t1 */
+    bnMulMod_(curve->t1, ptP->y, ptP->y, curve->p);       /* t1 = Y^2 */
+    bnMulMod_(curve->t0, ptP->x, mpiFour, curve->p);      /* t0 = 4 * X */
+    bnMulMod_(curve->S1, curve->t0, curve->t1, curve->p); /* S1 = t0 * t1 */
 
     /* M = 3*(X + Z^2)*(X - Z^2), use scratch variable U1 to store M value */
-    bnMulMod_(curve->t2, ptP->z, ptP->z, curve->p, curve);       /* t2 = Z^2 */
+    bnMulMod_(curve->t2, ptP->z, ptP->z, curve->p);       /* t2 = Z^2 */
     bnCopy(curve->t0, ptP->x);
-    bnAddMod_(curve->t0, curve->t2, curve->p);                   /* t0 = X + t2  */
-    bnMulMod_(curve->t3, curve->t0, mpiThree, curve->p, curve);  /* t3 = 3 * t0 */
+    bnAddMod_(curve->t0, curve->t2, curve->p);            /* t0 = X + t2  */
+    bnMulMod_(curve->t3, curve->t0, mpiThree, curve->p);  /* t3 = 3 * t0 */
     bnCopy(curve->t0, ptP->x);
-    bnSubMod_(curve->t0, curve->t2, curve->p);                   /* t0 = X - t2 */
-    bnMulMod_(curve->U1, curve->t3, curve->t0, curve->p, curve); /* M = t3 * t0 */
+    bnSubMod_(curve->t0, curve->t2, curve->p);            /* t0 = X - t2 */
+    bnMulMod_(curve->U1, curve->t3, curve->t0, curve->p); /* M = t3 * t0 */
     
     /* X' = M^2 - 2*S */
-    bnMulMod_(curve->t2, curve->U1, curve->U1, curve->p, curve); /* t2 = M^2 */
-    bnMulMod_(curve->t0, curve->S1, mpiTwo, curve->p, curve);    /* t0 = S * 2 */
+    bnMulMod_(curve->t2, curve->U1, curve->U1, curve->p); /* t2 = M^2 */
+    bnMulMod_(curve->t0, curve->S1, mpiTwo, curve->p);    /* t0 = S * 2 */
     bnCopy(R->x, curve->t2);
-    bnSubMod_(R->x, curve->t0, curve->p);                        /* X' = t2 - t0 */
+    bnSubMod_(R->x, curve->t0, curve->p);                 /* X' = t2 - t0 */
 
     /* Y' = M*(S - X') - 8*Y^4 */
-    bnMulMod_(curve->t3, curve->t1, curve->t1, curve->p, curve); /* t3 = Y^4 (t1 saved above) */
-    bnMulMod_(curve->t2, curve->t3, mpiEight, curve->p, curve);  /* t2 = t3 * 8 */
+    bnMulMod_(curve->t3, curve->t1, curve->t1, curve->p); /* t3 = Y^4 (t1 saved above) */
+    bnMulMod_(curve->t2, curve->t3, mpiEight, curve->p); /* t2 = t3 * 8 */
     bnCopy(curve->t3, curve->S1);
-    bnSubMod_(curve->t3, R->x, curve->p);                        /* t3 = S - X' */
-    bnMulMod_(curve->t0, curve->U1, curve->t3, curve->p, curve); /* t0 = M * t3 */
+    bnSubMod_(curve->t3, R->x, curve->p);                 /* t3 = S - X' */
+    bnMulMod_(curve->t0, curve->U1, curve->t3, curve->p); /* t0 = M * t3 */
     bnCopy(R->y, curve->t0);
-    bnSubMod_(R->y, curve->t2, curve->p);                        /* Y' = t0 - t2 */
+    bnSubMod_(R->y, curve->t2, curve->p);                 /* Y' = t0 - t2 */
 
     /* Z' = 2*Y*Z */
-    bnMulMod_(curve->t0, ptP->y, mpiTwo, curve->p, curve);       /* t0 = 2 * Y */
-    bnMulMod_(R->z, curve->t0, ptP->z, curve->p, curve);         /* Z' = to * Z */
+    bnMulMod_(curve->t0, ptP->y, mpiTwo, curve->p);       /* t0 = 2 * Y */
+    bnMulMod_(R->z, curve->t0, ptP->z, curve->p);         /* Z' = to * Z */
 
     if (P == R)
         FREE_EC_POINT(&tP);
@@ -622,72 +371,8 @@
     return ret;
 }
 
-static int ecDoublePointEd(const EcCurve *curve, EcPoint *R, const EcPoint *P)
-{
-    EcPoint tP;
-    const EcPoint *ptP = 0;
-
-    /* Check for overlapping arguments, copy if necessary and set pointer */
-    if (P == R) {
-        INIT_EC_POINT(&tP);
-        ptP = &tP;
-        bnCopy(tP.x, P->x);
-        bnCopy(tP.y, P->y);
-        bnCopy(tP.z, P->z);
-    }
-    else 
-        ptP = P;
-
-    /* Compute B, C, D, H, E */
-    bnCopy(curve->t1, ptP->x);
-    bnAddMod_(curve->t1, ptP->y, curve->p);
-    bnSquareMod_(curve->t0, curve->t1, curve->p, curve);       /* t0 -> B */
-
-    bnSquareMod_(R->x, ptP->x, curve->p, curve);               /* Rx -> C */
-
-    bnSquareMod_(R->y, ptP->y, curve->p, curve);               /* Ry -> D */
-
-    bnSquareMod_(R->z, ptP->z, curve->p, curve);               /* Rz -> H */
-    bnAddMod_(R->z, R->z, curve->p);                           /* Rz -> 2H */
-
-    bnCopy(curve->t1, R->x);
-    bnAddMod_(curve->t1, R->y, curve->p);                      /* t1 -> E */
-
-    /* Compute Ry */
-    bnCopy(curve->t2, R->x);
-    bnSubMod_(curve->t2, R->y, curve->p);                      /* C - D */
-    bnMulMod_(R->y, curve->t1, curve->t2, curve->p, curve);    /* E * t3; Ry */
-
-    /* Compute Rx */
-    bnSubMod_(curve->t0, curve->t1, curve->p);                 /* B - E; sub result */
-    bnCopy(curve->t2, curve->t1);
-    bnSubMod_(curve->t2, R->z, curve->p);                      /* t2 -> J; (E - 2H) */
-    bnMulMod_(R->x, curve->t2, curve->t0, curve->p, curve);    /* J * t0 */
-
-    /* Compute Rz */
-    bnMulMod_(R->z, curve->t2, curve->t1, curve->p, curve);    /* J * E */
-
-    if (P == R)
-        FREE_EC_POINT(&tP);
-
-    return 0;
-}
-
-/* 
- * Curve25519 has no specific Double Point function, all inside curve25519-donna
- */
-static int ecDoublePoint25519(const EcCurve *curve, EcPoint *R, const EcPoint *P)
-{
-    return -2;
-}
-
 /* Add two elliptic curve points. Any of them may be the same object. */
-int ecAddPoint(const EcCurve *curve, EcPoint *R, const EcPoint *P, const EcPoint *Q)
-{
-    return curve->addOp(curve, R, P, Q);
-}
-
-static int ecAddPointNist(const EcCurve *curve, EcPoint *R, const EcPoint *P, const EcPoint *Q)
+int ecAddPoint(const NistECpCurve *curve, EcPoint *R, const EcPoint *P, const EcPoint *Q)
 {
     int ret = 0;
 
@@ -739,23 +424,23 @@
         ptQ = Q;
 
     /* U1 = X1*Z2^2, where X1: P->x, Z2: Q->z */
-    bnMulMod_(curve->t1, ptQ->z, ptQ->z, curve->p, curve);    /* t1 = Z2^2 */
-    bnMulMod_(curve->U1, ptP->x, curve->t1, curve->p, curve); /* U1 = X1 * z_2 */
+    bnMulMod_(curve->t1, ptQ->z, ptQ->z, curve->p);    /* t1 = Z2^2 */
+    bnMulMod_(curve->U1, ptP->x, curve->t1, curve->p); /* U1 = X1 * z_2 */
 
     /* S1 = Y1*Z2^3, where Y1: P->y */
-    bnMulMod_(curve->t1, curve->t1, ptQ->z, curve->p, curve); /* t1 = Z2^3 */
-    bnMulMod_(curve->S1, ptP->y, curve->t1, curve->p, curve); /* S1 = Y1 * z_2 */
+    bnMulMod_(curve->t1, curve->t1, ptQ->z, curve->p); /* t1 = Z2^3 */
+    bnMulMod_(curve->S1, ptP->y, curve->t1, curve->p); /* S1 = Y1 * z_2 */
 
     /* U2 = X2*Z1^2, where X2: Q->x, Z1: P->z */
-    bnMulMod_(curve->t1, ptP->z, ptP->z, curve->p, curve);    /* t1 = Z1^2 */
-    bnMulMod_(curve->H, ptQ->x, curve->t1, curve->p, curve);  /* H = X2 * t1 (store U2 in H) */
+    bnMulMod_(curve->t1, ptP->z, ptP->z, curve->p);    /* t1 = Z1^2 */
+    bnMulMod_(curve->H, ptQ->x, curve->t1, curve->p);  /* H = X2 * t1 (store U2 in H) */
 
     /* H = U2 - U1 */
     bnSubMod_(curve->H, curve->U1, curve->p);
 
     /* S2 = Y2*Z1^3, where Y2: Q->y */
-    bnMulMod_(curve->t1, curve->t1, ptP->z, curve->p, curve); /* t1 = Z1^3 */
-    bnMulMod_(curve->R, ptQ->y, curve->t1, curve->p, curve);  /* R = Y2 * t1 (store S2 in R) */
+    bnMulMod_(curve->t1, curve->t1, ptP->z, curve->p); /* t1 = Z1^3 */
+    bnMulMod_(curve->R, ptQ->y, curve->t1, curve->p);  /* R = Y2 * t1 (store S2 in R) */
 
     /* R = S2 - S1 */
     bnSubMod_(curve->R, curve->S1, curve->p);
@@ -773,26 +458,26 @@
         return ecDoublePoint(curve, R, P);
     }
     /* X3 = R^2 - H^3 - 2*U1*H^2, where X3: R->x */
-    bnMulMod_(curve->t0, curve->H, curve->H, curve->p, curve);   /* t0 = H^2 */
-    bnMulMod_(curve->t1, curve->U1, curve->t0, curve->p, curve); /* t1 = U1 * t0, (hold t1) */
-    bnMulMod_(curve->t0, curve->t0, curve->H, curve->p, curve);  /* t0 = H^3, (hold t0) */
-    bnMulMod_(curve->t2, curve->R, curve->R, curve->p, curve);   /* t2 = R^2 */
+    bnMulMod_(curve->t0, curve->H, curve->H, curve->p);   /* t0 = H^2 */
+    bnMulMod_(curve->t1, curve->U1, curve->t0, curve->p); /* t1 = U1 * t0, (hold t1) */
+    bnMulMod_(curve->t0, curve->t0, curve->H, curve->p);  /* t0 = H^3, (hold t0) */
+    bnMulMod_(curve->t2, curve->R, curve->R, curve->p);   /* t2 = R^2 */
     bnCopy(curve->t3, curve->t2);
-    bnSubMod_(curve->t3, curve->t0, curve->p);                   /* t3 = t2 - t0, (-H^3)*/
-    bnMulMod_(curve->t2, mpiTwo, curve->t1, curve->p, curve);    /* t2 = 2 * t1 */
+    bnSubMod_(curve->t3, curve->t0, curve->p);            /* t3 = t2 - t0, (-H^3)*/
+    bnMulMod_(curve->t2, mpiTwo, curve->t1, curve->p);    /* t2 = 2 * t1 */
     bnCopy(R->x, curve->t3);
-    bnSubMod_(R->x, curve->t2, curve->p);                        /* X3 = t3 - t2 */
+    bnSubMod_(R->x, curve->t2, curve->p);                 /* X3 = t3 - t2 */
 
     /* Y3 = R*(U1*H^2 - X3) - S1*H^3, where Y3: R->y */
-    bnSubMod_(curve->t1, R->x, curve->p);                        /* t1 = t1 - X3, overwrites t1 now */
-    bnMulMod_(curve->t2, curve->R, curve->t1, curve->p, curve);  /* t2 = R * z_2 */
-    bnMulMod_(curve->S1, curve->S1, curve->t0, curve->p, curve); /* S1 = S1 * t0, (t0 has H^3) */
+    bnSubMod_(curve->t1, R->x, curve->p);                 /* t1 = t1 - X3, overwrites t1 now */
+    bnMulMod_(curve->t2, curve->R, curve->t1, curve->p);  /* t2 = R * z_2 */
+    bnMulMod_(curve->S1, curve->S1, curve->t0, curve->p); /* S1 = S1 * t0, (t0 has H^3) */
     bnCopy(R->y, curve->t2);
-    bnSubMod_(R->y, curve->S1, curve->p);                        /* Y3 = t2 - S1 */
+    bnSubMod_(R->y, curve->S1, curve->p);                 /* Y3 = t2 - S1 */
 
     /* Z3 = H*Z1*Z2, where Z1: P->z, Z2: Q->z, Z3: R->z */
-    bnMulMod_(curve->t2, curve->H, P->z, curve->p, curve);       /* t2 = H * Z1 */
-    bnMulMod_(R->z, curve->t2, Q->z, curve->p, curve);           /* Z3 = t2 * Z2 */
+    bnMulMod_(curve->t2, curve->H, P->z, curve->p);       /* t2 = H * Z1 */
+    bnMulMod_(R->z, curve->t2, Q->z, curve->p);           /* Z3 = t2 * Z2 */
 
     if (P == R)
         FREE_EC_POINT(&tP);
@@ -801,117 +486,10 @@
     return ret;
 }
 
-/*
- * Refer to the document: Faster addition and doubling on elliptic curves; Daniel J. Bernstein and Tanja Lange
- * section 4.
- *
- * This function is a variant of the 'addition'. The function returns the result in an own curve point
- * and does not overwrite its input parameters.
- */
-static int ecAddPointEd(const EcCurve *curve, EcPoint *R, const EcPoint *P, const EcPoint *Q)
+int ecMulPointScalar(const NistECpCurve *curve, EcPoint *R, const EcPoint *P, const BigNum *scalar)
 {
-    EcPoint tP, tQ;
-    const EcPoint *ptP = 0;
-    const EcPoint *ptQ = 0;
 
-    /* if P is (@,@), R = Q */
-    if (!bnCmp(P->z, mpiZero)) {
-        bnCopy(R->x, Q->x);
-        bnCopy(R->y, Q->y);
-        bnCopy(R->z, Q->z);
-        return 0;
-    }
-
-    /* if Q is (@,@), R = P */
-    if (!bnCmp(Q->z, mpiZero)) {
-        bnCopy(R->x, P->x);
-        bnCopy(R->y, P->y);
-        bnCopy(R->z, P->z);
-        return 0;
-    }
-
-    /* Check for overlapping arguments, copy if necessary and set pointers */
-    if (P == R) {
-        INIT_EC_POINT(&tP);
-        ptP = &tP;
-        bnCopy(tP.x, P->x);
-        bnCopy(tP.y, P->y);
-        bnCopy(tP.z, P->z);
-    }
-    else 
-        ptP = P;
-
-    if (Q == R) {
-        INIT_EC_POINT(&tQ);
-        ptQ = &tQ;
-        bnCopy(tQ.x, Q->x);
-        bnCopy(tQ.y, Q->y);
-        bnCopy(tQ.z, Q->z);
-    }
-    else
-        ptQ = Q;
-
-    /* Compute A, C, D first */
-    bnMulMod_(R->z, ptP->z, ptQ->z, curve->p, curve);            /* Rz -> A; (Z1 * z2); Rz becomes R3 */
-    bnMulMod_(R->x, ptP->x, ptQ->x, curve->p, curve);            /* Rx -> C; (X1 * X2); Rx becomes R1 */
-    bnMulMod_(R->y, ptP->y, ptQ->y, curve->p, curve);            /* Ry -> D; (Y1 * Y2); Ry becomes R2 */
-
-    /* Compute large parts of X3 equation, sub result in t0 */
-    bnCopy(curve->t0, ptP->x);
-    bnAddMod_(curve->t0, ptP->y, curve->p);                      /* t0 -> X1 + Y1 */
-    bnCopy(curve->t1, ptQ->x);
-    bnAddMod_(curve->t1, ptQ->y, curve->p);                      /* t1 -> X2 + Y2 */
-    bnMulMod_(curve->t2, curve->t0, curve->t1, curve->p, curve); /* t2 = t0 * t1 */
-    bnSubMod_(curve->t2, R->x, curve->p);                        /* t2 - C */
-    bnSubMod_(curve->t2, R->y, curve->p);                        /* t2 - D */
-    bnMulMod_(curve->t0, curve->t2, R->z, curve->p, curve);      /* t0 -> R7; (t2 * A); sub result */
-
-    /* Compute E */
-    bnMulMod_(curve->t2, R->x, R->y, curve->p, curve);           /* t2 = C * D */
-    bnMulMod_(curve->t1, curve->t2, curve->a, curve->p, curve);  /* t1 -> E; t1 new R8 */
-
-    /* Compute part of Y3 equation, sub result in t2 */
-    bnSubMod_(R->y, R->x, curve->p);                             /* Ry = D - C; sub result */
-    bnMulMod_(curve->t2, R->y, R->z, curve->p, curve);           /* t2 = Ry * A; sub result */
-
-    /* Compute B */
-    bnSquareMod_(R->z, R->z, curve->p, curve);                   /* Rz -> B; (A^2) */
-
-    /* Compute F */
-    bnCopy(curve->t3, R->z);
-    bnSubMod_(curve->t3, curve->t1, curve->p);                   /* t3 -> F; (B - E) */
-
-    /* Compute G */
-    bnAddMod_(R->z, curve->t1, curve->p);                        /* Rz -> G; (B + E) */
-
-    /* Compute, X, Y, Z results */
-    bnMulMod_(R->x, curve->t3, curve->t0, curve->p, curve);      /* Rx = F * t0 */
-    bnMulMod_(R->y, curve->t2, R->z, curve->p, curve);           /* Ry = t2 * G */
-    bnMulMod_(R->z, curve->t3, R->z, curve->p, curve);           /* Rz = F * G */
-
-    if (P == R)
-        FREE_EC_POINT(&tP);
-    if (Q == R)
-        FREE_EC_POINT(&tQ);
-
-    return 0;
-}
-
-/* 
- * Curve25519 has no specific Add Point function, all inside curve25519-donna
- */
-static int ecAddPoint25519(const EcCurve *curve, EcPoint *R, const EcPoint *P, const EcPoint *Q)
-{
-    return -2;
-}
-
-int ecMulPointScalar(const EcCurve *curve, EcPoint *R, const EcPoint *P, const BigNum *scalar)
-{
-    return curve->mulScalar(curve, R, P, scalar);
-}
-
-static int ecMulPointScalarNormal(const EcCurve *curve, EcPoint *R, const EcPoint *P, const BigNum *scalar)
-{
+    /* MPI_CHK below macro requires a 'ret' variable and a cleanup label */
     int ret = 0;
     int i;
     int bits = bnBits(scalar);
@@ -937,26 +515,7 @@
     return ret;
 }
 
-/* 
- * This function uses BigNumber only as containers to transport the 32 byte data.
- * This makes it compliant to the other functions and thus higher-level API does not change.
- * 
- * curve25519_donna function uses data in little endian format.
- */
-static int ecMulPointScalar25519(const EcCurve *curve, EcPoint *R, const EcPoint *P, const BigNum *scalar)
-{
-    uint8_t basepoint[32], secret[32], result[32];
-
-    bnExtractLittleBytes(P->x, basepoint, 0, 32);  /* 25519 function requires the X coordinate only (compressed) */
-    bnExtractLittleBytes(scalar, secret, 0, 32);
-    curve25519_donna(result, secret, basepoint);
-    bnInsertLittleBytes(R->x, result, 0, 32);
-    return 0;
-}
-
 #ifdef WEAKRANDOM
-#include <fcntl.h>
-
 /*
  * A standard random number generator that uses the portable random() system function.
  *
@@ -964,15 +523,10 @@
  */
 static int _random(unsigned char *output, size_t len)
 {
-    size_t num = 0;
+    size_t i;
 
-    int rnd = open("/dev/urandom", O_RDONLY);
-    if (rnd >= 0) {
-        num = read(rnd, output, len);
-        close(rnd);
-    }
-    else
-        return num;
+    for(i = 0; i < len; ++i )
+        output[i] = random();
 
     return( 0 );
 }
@@ -984,12 +538,7 @@
 }
 #endif
 
-int ecGenerateRandomNumber(const EcCurve *curve, BigNum *d)
-{
-    return curve->randomOp(curve, d);
-}
-
-static int ecGenerateRandomNumberNist(const EcCurve *curve, BigNum *d)
+int ecGenerateRandomNumber(const NistECpCurve *curve, BigNum *d)
 {
     BigNum c, nMinusOne;
 
@@ -1019,677 +568,3 @@
 
     return 0;
 }
-
-static int ecGenerateRandomNumber3617(const EcCurve *curve, BigNum *d)
-{
-    unsigned char random[52];
-    _random(random, 52);
-
-    /* prepare the secret random data: clear bottom 3 bits. Clearing top 2 bits
-     * makes is a 414 bit value
-     */
-    random[51] &= ~0x7;
-    random[0] &= 0x3f;
-    /* convert the random data into big numbers */
-    bnInsertBigBytes(d, random, 0, 52);
-    return 0;
-}
-
-static int ecGenerateRandomNumber25519(const EcCurve *curve, BigNum *d)
-{
-    unsigned char random[32];
-    _random(random, 32);
-
-    /* No specific preparation. The curve25519_donna functions prepares the data.
-     *
-     * convert the random data into big numbers. the bigNumber is a container only.
-     * we don not use the big number for any arithmetic
-     */
-    bnInsertLittleBytes(d, random, 0, 32);
-    return 0;
-
-}
-
-int ecCheckPubKey(const EcCurve *curve, const EcPoint *pub)
-{
-    return curve->checkPubOp(curve, pub);
-}
-
-static int ecCheckPubKeyNist(const NistECpCurve *curve, const EcPoint *pub)
-{
-    /* Represent point at infinity by (0, 0), make sure it's not that */
-    if (bnCmpQ(pub->x, 0) == 0 && bnCmpQ(pub->y, 0) == 0) {
-        return 0;
-    }
-    /* Check that coordinates are within range */
-    if (bnCmpQ(pub->x, 0) < 0 || bnCmp(pub->x, curve->p) >= 0) {
-        return 0;
-    }
-    if (bnCmpQ(pub->y, 0) < 0 || bnCmp(pub->y, curve->p) >= 0) {
-        return 0;
-    }
-    /* Check that point satisfies EC equation y^2 = x^3 - 3x + b, mod P */
-    bnSquareMod_(curve->t1, pub->y, curve->p, curve);
-    bnSquareMod_(curve->t2, pub->x, curve->p, curve);
-    bnSubQMod_(curve->t2, 3, curve->p);
-    bnMulMod_(curve->t2, curve->t2, pub->x, curve->p, curve);
-    bnAddMod_(curve->t2, curve->b, curve->p);
-    if (bnCmp (curve->t1, curve->t2) != 0) {
-        return 0;
-    }
-    return 1;
-
-}
-
-static int ecCheckPubKey3617(const EcCurve *curve, const EcPoint *pub)
-{
-    /* Represent point at infinity by (0, 0), make sure it's not that */
-    if (bnCmpQ(pub->x, 0) == 0 && bnCmpQ(pub->y, 0) == 0) {
-        return 0;
-    }
-    /* Check that coordinates are within range */
-    if (bnCmpQ(pub->x, 0) < 0 || bnCmp(pub->x, curve->p) >= 0) {
-        return 0;
-    }
-    if (bnCmpQ(pub->y, 0) < 0 || bnCmp(pub->y, curve->p) >= 0) {
-        return 0;
-    }
-    /* Check that point satisfies EC equation x^2+y^2 = 1+3617x^2y^2, mod P */
-    bnSquareMod_(curve->t1, pub->y, curve->p, curve);
-    bnSquareMod_(curve->t2, pub->x, curve->p, curve);
-    bnCopy(curve->t3, curve->t1);                                /* Load t3 */
-    bnAddMod_(curve->t3, curve->t2, curve->p);                   /* t3 = t1 + t2, (x^2+y^2)*/
-
-    bnMulMod_(curve->t0, curve->a, curve->t1, curve->p, curve);  /* t0 = a * t1,  (3617 * x^2) */
-    bnMulMod_(curve->t0, curve->t0, curve->t2, curve->p, curve); /* t0 = t0 * t1, (3617 * x^2 * y^2) */
-    bnAddMod_(curve->t0, mpiOne, curve->p);                      /* t0 = t0 + 1,  (3617 * x^2 * y^2 + 1) */
-
-    if (bnCmp (curve->t0, curve->t3) != 0) {
-        return 0;
-    }
-    return 1;
-}
-
-/**
- * According to http://cr.yp.to/ecdh.html#validate no validation is required if used for Diffie-Hellman
- * thus always return success.
- */
-static int ecCheckPubKey25519(const EcCurve *curve, const EcPoint *pub)
-{
-    return 1;
-}
-
-static int mod3617(BigNum *r, const BigNum *a, const BigNum *modulo)
-{
-    unsigned char buffer[52] = {0};
-    int cmp;
-    BigNum tmp;
-
-    bnBegin(&tmp);
-    cmp = bnCmp(modulo, a);
-    if (cmp == 0) {             /* a is equal modulo, set resul to zero */
-        bnSetQ(r, 0);
-        return 0;
-    }
-    if (cmp > 0) {              /* modulo is greater than a - copy a to r and return it */
-        bnCopy(r, a);
-        return 0;
-    }
-    bnExtractLittleBytes(a, buffer, 0, 52);
-    buffer[51] &= 0x3f;
-
-    bnCopy(&tmp, a);
-    bnRShift(&tmp, 414);
-    bnCopy(r, &tmp);
-    bnLShift(r, 4);
-    bnAdd(r, &tmp);
-
-    bnInsertLittleBytes(&tmp, buffer, 0, 52);
-
-    bnAdd(r, &tmp);
-    while (bnCmp(r, modulo) >= 0) {
-        bnSub(r, modulo);
-    }
-    bnEnd(&tmp);
-    return 0;
-}
-
-/* 
- * Curve25519 has no specific modulo function, all inside curve25519-donna
- */
-static int mod25519(BigNum *r, const BigNum *a, const BigNum *modulo)
-{
-    return -2;
-}
-
-/*
- * Beware: Here are the dragons.
- *
- * The modulo implementations for the NIST curves. For more detailled information see
- * FIPS 186-3, chapter D.2 and other papers about Generailzed Mersenne numbers.
- *
- * I use byte operations to perfom the additions with carry. On a little endian machine
- * this saves conversion from/to big endian format if I would use integers for example. Also
- * using byte addition into a short carry accumulator works on every word size and avoids
- * complex testing and handling of wordsizes and big/little endian stuff.
- *
- */
-
-/* new modulo for 192bit curve */
-static int newMod192(BigNum *r, const BigNum *a, const BigNum *modulo)
-{
-    unsigned char buffer[200] = {0};
-    unsigned char *pt;
-    unsigned char *ps1;
-    unsigned char *ps2;
-    unsigned char *ps3;
-    short ac;
-    int cmp;
-
-    /* Binary big number representation in PolarSSL is always big endian
-     *
-     * the least significant 64bit large word starts at byte offset 40,
-     * the least significant 32bit word starts at byte offset 44
-     * the least significant byte starts at byte offset 47
-     *
-     *           S3    S2   S1          T
-     *                            /-----^------\
-     *           A5    A4   A3    A2    A1    A0
-     * 64bit  0     1     2     3     4     5
-     *        |--+--|--+--|--+--|--+--|--+--|--+--|
-     * 32bit  0  1  2  3  4  5  6  7  8  9 10 11
-     *
-     * perform T + S1 + S2 + S3 mod p
-
-     * where T  = (A2 || A1 || A0)
-     *     + S1 = ( 0 || A3 || A3)
-     *     + S2 = (A4 || A4 ||  0)
-     *     + S3 = (A5 || A5 || A5)
-     *
-     * TODO: error check if input variable is > modulo^2 (do normal mpi_mod_mpi),
-     */
-
-    /* TODO: check if a is > modulo^2 */
-    cmp = bnCmp(modulo, a);
-    if (cmp == 0) {             /* a is equal modulo, set resul to zero */
-        bnSetQ(r, 0);
-        return 0;
-    }
-    if (cmp > 0) {              /* modulo is greater than a - copy a to r and return it */
-        bnCopy(r, a);
-        return 0;
-    }
-    bnExtractBigBytes(a, buffer, 0, bnBytes(modulo)*2);
-
-    /* 6 'A' words, each word is 8 byte. Compute offset to least significant byte of word X */
-#define A(X) buffer + (((6-X)*8)-1)
-
-    ac = 0;
-
-    pt = A(0);      /* pt points to least significant byte of A0  */
-
-    /* Add up first 8 byte word, no need to add ps2 */
-    ps1 = A(3);        /* ps1 points to least significant byte of S1 (A3) */
-    ps3 = A(5);        /* ps3 points to least significant byte of S3 (A5)*/
-
-    /* Each block processes one 32 bit word, big endian, using byte operations */
-    ac += *pt + *ps1--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-
-    ac += *pt + *ps1--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-
-    /* Add up second 8 byte word, all three S words are used here */
-    ps1 = A(3); ps2 = A(4); ps3 = A(5);
-
-    ac += *pt + *ps1--; ac += *ps2--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1--; ac += *ps2--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1--; ac += *ps2--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1--; ac += *ps2--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-
-    ac += *pt + *ps1--; ac += *ps2--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1--; ac += *ps2--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1--; ac += *ps2--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1--; ac += *ps2--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-
-    /* Add up third 8 byte word, no need to add S1 word */
-    ps2 = A(4); ps3 = A(5);
-
-    ac += *pt + *ps2--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-
-    ac += *pt + *ps2--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; *pt-- = ac; ac >>= 8;
-
-    /* In this function we cannot have a negative carry and at most a carry of 2
-     * thus just subtract the modulo until we are less than modulo
-     */
-    bnSetQ(r, 0);
-
-    *(A(3)) = ac;      /* Store the carry */
-    bnInsertBigBytes(r, A(3), 0, 25);  /* 25: 3 * 8 byte words + 1 carry byte */
-    while (bnCmp(r, modulo) >= 0) {
-        bnSub(r, modulo);
-    }
-    return 0;
-}
-#undef A
-
-/* new modulo for 256bit curve */
-static int newMod256(BigNum *r, const BigNum *a, const BigNum *modulo)
-{
-    unsigned char buffer[200] = {0};
-    unsigned char *pt;
-    unsigned char *ps1;
-    unsigned char *ps2;
-    unsigned char *ps3;
-    unsigned char *ps4;
-
-    unsigned char *pd1;
-    unsigned char *pd2;
-    unsigned char *pd3;
-    unsigned char *pd4;
-    short ac;
-    int cmp;
-
-    /* Binary big number representation in PolarSSL is always big endian
-     *
-     * the least significant byte starts at byte offset 63
-     *
-     *                                                                    T
-     *                                                  /-----------------^------------------\
-     *          A15  A14  A13  A12  A11  A10  A9   A8   A7   A6   A5   A4   A3   A2   A1   A0
-     *        |----+----|----+----|----+----|----+----|----+----|----+----|----+----|----+----|
-     * offset 0    4    8   12   16   20   24   28   32   36   40   44   48   52   56   60    64
-     *
-     * T  = (  A7 ||  A6 ||  A5 ||  A4 ||  A3 ||  A2 ||  A1 ||  A0 )
-     *
-     * S1 = ( A15 || A14 || A13 || A12 || A11 ||  00 ||  00 ||  00 )
-     * S2 = (  00 || A15 || A14 || A13 || A12 ||  00 ||  00 ||  00 )
-     * S3 = ( A15 || A14 ||  00 ||  00 ||  00 || A10 ||  A9 ||  A8 )
-     * S4 = (  A8 || A13 || A15 || A14 || A13 || A11 || A10 ||  A9 )
-     * D1 = ( A10 ||  A8 ||  00 ||  00 ||  00 || A13 || A12 || A11 )
-     * D2 = ( A11 ||  A9 ||  00 ||  00 || A15 || A14 || A13 || A12 )
-     * D3 = ( A12 ||  00 || A10 ||  A9 ||  A8 || A15 || A14 || A13 )
-     * D4 = ( A13 ||  00 || A11 || A10 ||  A9 ||  00 || A15 || A14 )
-     *
-     * perform B = T + 2*S1 + 2*S2 + S3 + S4 - D1 - D2 - D3 - D4 mod p
-     *
-     * TODO: error check if input variable is > modulo^2 (do normal mpi_mod_mpi),
-     */
-
-    cmp = bnCmp(modulo, a);
-    if (cmp == 0) {             /* a is equal modulo, set resul to zero */
-        bnSetQ(r, 0);
-        return 0;
-    }
-    if (cmp > 0) {              /* modulo is greater than a - copya to r and return it */
-        bnCopy(r, a);
-        return 0;
-    }
-    bnExtractBigBytes(a, buffer, 0, bnBytes(modulo)*2);
-
-    /* 16 'A' words, each word is 4 byte. Compute offset to least significant byte of word X */
-#define A(X) buffer + (((16-X)*4)-1)
-
-    ac = 0;
-
-    pt = A(0);          /* pt points to least significant byte of A0  */
-
-    /* Set up to add up data that goes into A0 (right-most column abover); S1, S2 not used */
-    ps3 = A(8);         /* ps3 points to least significant byte of S3 */
-    ps4 = A(9);         /* ps4 points to least significant byte of S4 */
-    pd1 = A(11);        /* pd1 points to least significant byte of D1 */
-    pd2 = A(12);        /* pd2 points to least significant byte of D2 */
-    pd3 = A(13);        /* pd3 points to least significant byte of D3 */
-    pd4 = A(14);        /* pd4 points to least significant byte of D4 */
-
-    /* Each block processes one 32 bit word, big endian, using byte operations */
-    ac += *pt + *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-
-    /* Set up to add up data that goes into A1; S1, S2 not used */
-    ps3 = A(9);  ps4 = A(10); pd1 = A(12); pd2 = A(13); pd3 = A(14); pd4 = A(15);
-    ac += *pt + *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-
-    /* Set up to add up data that goes into A2; S1, S2, D4 not used */
-    ps3 = A(10); ps4 = A(11); pd1 = A(13); pd2 = A(14); pd3 = A(15);
-    ac += *pt + *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; ac -= *pd3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; ac -= *pd3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; ac -= *pd3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; ac -= *pd3--; *pt-- = ac; ac >>= 8;
-
-    /* Set up to add up data that goes into A3; S3, D1 not used */
-    ps1 = A(11); ps2 = A(12); ps4 = A(13); pd2 = A(15); pd3 = A(8); pd4 = A(9);
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps2;ac += *ps2--;  ac += *ps4--; ac -= *pd2--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps2;ac += *ps2--;  ac += *ps4--; ac -= *pd2--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps2;ac += *ps2--;  ac += *ps4--; ac -= *pd2--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps2;ac += *ps2--;  ac += *ps4--; ac -= *pd2--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-
-    /* Set up to add up data that goes into A4; S3, D1, D2 not used */
-    ps1 = A(12); ps2 = A(13); ps4 = A(14); pd3 = A(9); pd4 = A(10);
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps2;ac += *ps2--;  ac += *ps4--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps2;ac += *ps2--;  ac += *ps4--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps2;ac += *ps2--;  ac += *ps4--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps2;ac += *ps2--;  ac += *ps4--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-
-    /* Set up to add up data that goes into A5; S3, D1, D2 not used */
-    ps1 = A(13); ps2 = A(14); ps4 = A(15); pd3 = A(10); pd4 = A(11);
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps2;ac += *ps2--;  ac += *ps4--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps2;ac += *ps2--;  ac += *ps4--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps2;ac += *ps2--;  ac += *ps4--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps2;ac += *ps2--;  ac += *ps4--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-
-    /* Set up to add up data that goes into A6; D3, D4 not used */
-    ps1 = A(14); ps2 = A(15); ps3 = A(14); ps4 = A(13); pd1 = A(8); pd2 = A(9);
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps2;ac += *ps2--;  ac += *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps2;ac += *ps2--;  ac += *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps2;ac += *ps2--;  ac += *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps2;ac += *ps2--;  ac += *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; *pt-- = ac; ac >>= 8;
-
-    /* Set up to add up data that goes into A7; S2 not used */
-    ps1 = A(15); ps3 = A(15); ps4 = A(8); pd1 = A(10); pd2 = A(11); pd3 = A(12); pd4 = A(13);
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--;  ac += *ps3--; ac += *ps4--; ac -= *pd1--; ac -= *pd2--; ac -= *pd3--; ac -= *pd4--; *pt-- = ac; ac >>= 8;
-
-    bnSetQ(r, 0);
-    if (ac > 0) {
-        *(A(8)) = ac;      /* Store the carry */
-        bnInsertBigBytes(r, A(8), 0, 33);  /* 33: 8 * 4 byte words + 1 carry byte */
-    }
-    /* Negative carry requires that we add the modulo (carry * -1) times to make
-     * the result positive. Then get the result mod(256).
-     */
-    else if (ac < 0) {
-        int msb, maxMsb;
-
-        *(A(8)) = 0;
-        bnInsertBigBytes(r, A(8), 0, 33);  /* 33: 8 * 4 byte words + 1 carry byte */
-        ac *= -1;
-        while (ac--) {
-            bnAdd(r, modulo);
-        }
-        maxMsb =  bnBits(modulo);
-        msb = bnBits(r) - maxMsb;
-        /* clear all bits above bit length of modulo. This length is 256 here, thus
-         * we effectiviely doing a mod(256)
-         */
-        if (msb > 0) {
-            BigNum tmp;
-            bnBegin(&tmp);
-            bnSetQ (&tmp, 1);
-            bnLShift (&tmp, maxMsb);
-            bnMod(r, r, &tmp);
-            bnEnd(&tmp);
-        }
-    }
-    else {
-        *(A(8)) = 0;
-        bnInsertBigBytes(r, A(8), 0, 33);  /* 33: 8 * 4 byte words + 1 carry byte */
-    }
-    while (bnCmp(r, modulo) >= 0) {
-        bnSub(r, modulo);
-    }
-    return 0;
-}
-#undef A
-
-
-/* new modulo for 384bit curve */
-static int newMod384(BigNum *r, const BigNum *a, const BigNum *modulo)
-{
-    unsigned char buffer[200] = {0};
-    unsigned char *pt;
-    unsigned char *ps1;
-    unsigned char *ps2;
-    unsigned char *ps3;
-    unsigned char *ps4;
-    unsigned char *ps5;
-    unsigned char *ps6;
-
-    unsigned char *pd1;
-    unsigned char *pd2;
-    unsigned char *pd3;
-    short ac;
-    int cmp;
-
-    /*
-     *
-     * the least significant byte starts at byte offset 97
-     *
-     *                                                                    T
-     *                                        /---------------------------^----------------------------\
-     *      A23 ......... A15  A14  A13  A12  A11  A10  A9   A8   A7   A6   A5   A4   A3   A2   A1   A0
-     *    |----+ ...... |----+----|----+----|----+----|----+----|----+----|----+----|----+----|----+----|
-     *
-     * T  = (A11 || A10 ||  A9 ||  A8 ||  A7 ||  A6 ||  A5 ||  A4 ||  A3 ||  A2 ||  A1 ||  A0)
-
-     * S1 = ( 00 ||  00 ||  00 ||  00 ||  00 || A23 || A22 || A21 ||  00 ||  00 ||  00 ||  00)
-     * S2 = (A23 || A22 || A21 || A20 || A19 || A18 || A17 || A16 || A15 || A14 || A13 || A12)
-     * S3 = (A20 || A19 || A18 || A17 || A16 || A15 || A14 || A13 || A12 || A23 || A22 || A21)
-     * S4 = (A19 || A18 || A17 || A16 || A15 || A14 || A13 || A12 || A20 ||  00 || A23 ||  00)
-     * S5 = ( 00 ||  00 ||  00 ||  00 || A23 || A22 || A21 || A20 ||  00 ||  00 ||  00 ||  00)
-     * S6 = ( 00 ||  00 ||  00 ||  00 ||  00 ||  00 || A23 || A22 || A21 ||  00 ||  00 || A20)
-     * D1 = (A22 || A21 || A20 || A19 || A18 || A17 || A16 || A15 || A14 || A13 || A12 || A23)
-     * D2 = ( 00 ||  00 ||  00 ||  00 ||  00 ||  00 ||  00 || A23 || A22 || A21 || A20 ||  00)
-     * D3 = ( 00 ||  00 ||  00 ||  00 ||  00 ||  00 ||  00 || A23 || A23 ||  00 ||  00 ||  00)
-     *
-     * perform B =  T + 2S1 + S2 + S3 + S4 + S5 + S6 – D1 – D2 – D3 mod p
-     *
-     * TODO: error check if input variable is > modulo^2 (do normal mpi_mod_mpi),
-     *       optimize if input is already < modulo  (just copy over in this case).
-     */
-
-    cmp = bnCmp(modulo, a);
-    if (cmp == 0) {             /* a is equal modulo, set resul to zero */
-        bnSetQ(r, 0);
-        return 0;
-    }
-    if (cmp > 0) {              /* modulo is greater than a - copy a to r and return it */
-        bnCopy(r, a);
-        return 0;
-    }
-
-    bnExtractBigBytes(a, buffer, 0, bnBytes(modulo)*2);
-
-    /* 24 'A' words, each word is 4 byte. Compute offset to least significant byte of word X */
-#define A(X) buffer + (((24-X)*4)-1)
-
-    ac = 0;
-
-    pt = A(0);      /* pt points to least significant byte of A0  */
-
-    /* Set up to add data that goes into A0; S1, S4, S5, D2, D3 not used */
-    ps2 = A(12); ps3 = A(21); ps6 = A(20); pd1 = A(23);
-
-    /* Each block processes one 32 bit word, big endian, using byte operations */
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps6--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps6--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps6--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps6--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-
-    /* Set up to add data that goes into A1; S1, S5, S6, D3 not used */
-    ps2 = A(13); ps3 = A(22); ps4 = A(23); pd1= A(12); pd2 = A(20);
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--;  ac -= *pd2--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--;  ac -= *pd2--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--;  ac -= *pd2--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--;  ac -= *pd2--; *pt-- = ac; ac >>= 8;
-
-    /* Set up to add data that goes into A2; S1, S4, S5, S6, D3 not used */
-    ps2 = A(14); ps3 = A(23); pd1 = A(13); pd2 = A(21);
-    ac += *pt + *ps2--; ac += *ps3--; ac -= *pd1--;  ac -= *pd2--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac -= *pd1--;  ac -= *pd2--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac -= *pd1--;  ac -= *pd2--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac -= *pd1--;  ac -= *pd2--; *pt-- = ac; ac >>= 8;
-
-    /* Set up to add data that goes into A3; S1, S5, S6 not used */
-    ps2 = A(15); ps3 = A(12); ps4 = A(20); ps6 = A(21); pd1 = A(14); pd2 = A(22); pd3 = A(23);
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps6--; ac -= *pd1--;  ac -= *pd2--; ac -= *pd3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps6--; ac -= *pd1--;  ac -= *pd2--; ac -= *pd3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps6--; ac -= *pd1--;  ac -= *pd2--; ac -= *pd3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps6--; ac -= *pd1--;  ac -= *pd2--; ac -= *pd3--; *pt-- = ac; ac >>= 8;
-
-    /* Set up to add data that goes into A4 */
-    ps1 = A(21); ps2 = A(16); ps3 = A(13); ps4 = A(12); ps5 = A(20); ps6 = A(22); pd1 = A(15); pd2 = A(23), pd3 = A(23);
-    ac += *pt + *ps1;ac += *ps1--; ac += *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps5--; ac += *ps6--; ac -= *pd1--;  ac -= *pd2--; ac -= *pd3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--; ac += *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps5--; ac += *ps6--; ac -= *pd1--;  ac -= *pd2--; ac -= *pd3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--; ac += *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps5--; ac += *ps6--; ac -= *pd1--;  ac -= *pd2--; ac -= *pd3--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--; ac += *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps5--; ac += *ps6--; ac -= *pd1--;  ac -= *pd2--; ac -= *pd3--; *pt-- = ac; ac >>= 8;
-
-    /* Set up to add data that goes into A5; D2, D3 not used */
-    ps1 = A(22); ps2 = A(17); ps3 = A(14); ps4 = A(13); ps5 = A(21); ps6 = A(23); pd1 = A(16);
-    ac += *pt + *ps1;ac += *ps1--; ac += *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps5--; ac += *ps6--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--; ac += *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps5--; ac += *ps6--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--; ac += *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps5--; ac += *ps6--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--; ac += *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps5--; ac += *ps6--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-
-    /* Set up to add data that goes into A6; S6, D2, D3 not used */
-    ps1 = A(23); ps2 = A(18); ps3 = A(15); ps4 = A(14); ps5 = A(22); pd1 = A(17);
-    ac += *pt + *ps1;ac += *ps1--; ac += *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps5--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--; ac += *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps5--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--; ac += *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps5--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps1;ac += *ps1--; ac += *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps5--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-
-    /* Set up to add data that goes into A7; S1, S6, D2, D3 not used */
-    ps2 = A(19); ps3 = A(16); ps4 = A(15); ps5 = A(23); pd1 = A(18);
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps5--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps5--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps5--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac += *ps5--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-
-    /* Set up to add data that goes into A8; S1, S5, S6, D2, D3 not used */
-    ps2 = A(20); ps3 = A(17); ps4 = A(16); pd1 = A(19);
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-
-    /* Set up to add data that goes into A9; S1, S5, S6, D2, D3 not used */
-    ps2 = A(21); ps3 = A(18); ps4 = A(17); pd1 = A(20);
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-
-    /* Set up to add data that goes into A10; S1, S5, S6, D2, D3 not used */
-    ps2 = A(22); ps3 = A(19); ps4 = A(18); pd1 = A(21);
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-
-    /* Set up to add data that goes into A10; S1, S5, S6, D2, D3 not used */
-    ps2 = A(23); ps3 = A(20); ps4 = A(19); pd1 = A(22);
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-    ac += *pt + *ps2--; ac += *ps3--; ac += *ps4--; ac -= *pd1--; *pt-- = ac; ac >>= 8;
-
-    bnSetQ(r, 0);
-    if (ac > 0) {
-        *(A(12)) = ac;      /* Store the carry */
-        bnInsertBigBytes(r, A(12), 0, 49);  /* 49: 12 * 4 byte words + 1 carry byte */
-    }
-    /* Negative carry requires that we add the modulo (carry * -1) times to make
-     * the result positive. Then get the result mod(256).
-     */
-    else if (ac < 0) {
-        int msb, maxMsb;
-
-        *(A(12)) = 0;
-        bnInsertBigBytes(r, A(12), 0, 49);  /* 49: 12 * 4 byte words + 1 carry byte */
-        ac *= -1;
-        while (ac--) {
-            bnAdd(r, modulo);
-        }
-        maxMsb =  bnBits(modulo);
-        msb = bnBits(r) - maxMsb;
-        /* clear all bits above bit length of modulo. This length is 384 here, thus
-         * we effectiviely doing a mod(384)
-         */
-        if (msb > 0) {
-            BigNum tmp;
-            bnBegin(&tmp);
-            bnSetQ (&tmp, 1);
-            bnLShift (&tmp, maxMsb);
-            bnMod(r, r, &tmp);
-            bnEnd(&tmp);
-        }
-    }
-    else {
-        *(A(12)) = 0;
-        bnInsertBigBytes(r, A(12), 0, 49);  /* 49: 12 * 4 byte words + 1 carry byte */
-    }
-    while (bnCmp(r, modulo) >= 0) {
-        bnSub(r, modulo);
-    }
-    return 0;
-}
-#undef A
-
-
-/* new modulo for 521bit curve, much easier because the prime for 521 is a real Mersenne prime */
-static int newMod521(BigNum *r, const BigNum *a, const BigNum *modulo)
-{
-    unsigned char buf1[200] = {0};
-    unsigned char buf2[200] = {0};
-    unsigned char *p1;
-    unsigned char *p2;
-    size_t modSize;
-    short ac = 0;
-    unsigned int i;
-    int cmp;
-
-    /* TODO: check if a is > modulo^2 */
-#if 0
-    if (a->s < 0)               /* is it a negative value? */
-        return bnMod(r, a, modulo);
-#endif
-    cmp = bnCmp(modulo, a);
-    if (cmp == 0) {             /* a is equal modulo, set resul to zero */
-        bnSetQ(r, 0);
-        return 0;
-    }
-    bnCopy(r, a);
-    if (cmp > 0) {              /* modulo is greater than a - return the prepared r */
-        return 0;
-    }
-    modSize = bnBytes(modulo);
-
-    bnExtractBigBytes(a, buf1, 0, modSize*2); /* a must be less modulo^2 */
-    buf1[modSize] &= 1;                   /* clear all bits except least significat */
-
-    bnRShift(r, 521);
-    bnExtractBigBytes(r, buf2, 0, modSize*2);
-    buf2[modSize] &= 1;
-
-    p1 = &buf2[131];            /* p1 is pointer to A0 */
-    p2 = &buf1[131];            /* p2 is pointer to A1 */
-
-    for (i = 0; i < modSize; i++) {
-        ac += *p1 + *p2--; *p1-- = ac; ac >>= 8;
-    }
-    bnSetQ(r, 0);
-    bnInsertBigBytes(r, p1+1, 0, modSize);
-
-    while (bnCmp(r, modulo) >= 0) {
-        bnSub(r, modulo);
-    }
-    return 0;
-}
-
diff --git a/jni/libzrtp/sources/bnlib/ec/ec.h b/jni/libzrtp/sources/bnlib/ec/ec.h
index 172ffd8..b1bb47c 100644
--- a/jni/libzrtp/sources/bnlib/ec/ec.h
+++ b/jni/libzrtp/sources/bnlib/ec/ec.h
@@ -29,34 +29,18 @@
     NIST224P = 2,
     NIST256P = 3,
     NIST384P = 4,
-    NIST521P = 5,
-    Curve25519 = 10,
-    Curve3617  = 11
-} Curves;
+    NIST521P = 5
+} NistCurves;
 
 /**
- * \brief This structure contains the x, y affine coordinates and the z value if we
- *        use projective coordinates during EC point arithmetic.
- */
-typedef struct _EcPoint {
-    BigNum *x, *y, *z;
-    BigNum tx, ty, tz;
-} EcPoint;
-
-/**
- * @brief This structure contains the value of EC curves over Prime Fields.
+ * @brief This structure contains the value of NIST EC curves over Prime Fields.
  *
- * The for NIST curves the field names correspond to the variable names defined in 
- * NIST FIPS 186-3, E.1.2. The <b>a</b> curve parameter is the constant -3 and is 
- * computed during initialization of the curve structure.
+ * The <b>a</b> curve parameter is the constant -3 and is computed during initialization
+ * of the curve structure.
  *
- * For other curves, for example curve3917 we have less parameters to fill in, mostly
- * the prime number, the base point, etc. Refer to the curve's initialization function
- * about the use of the fileds.
+ * The field names correspond to the variable names defined in NIST FIPS 186-3, E.1.2
  */
-struct EcCurve;
-struct EcCurve {
-    Curves id;
+typedef struct {
     BigNum _p;
     BigNum _n;
     BigNum _SEED;
@@ -78,18 +62,16 @@
        avoid to much memory allocation/deallocatio0n overhead */
   BigNum _S1, _U1, _H, _R, _t0, _t1, _t2, _t3;
   BigNum *S1, *U1, *H, *R, *t0, *t1, *t2, *t3;
-  int (*affineOp)(const struct EcCurve *curve, EcPoint *R, const EcPoint *P);
-  int (*doubleOp)(const struct EcCurve *curve, EcPoint *R, const EcPoint *P);
-  int (*addOp)(const struct EcCurve *curve, EcPoint *R, const EcPoint *P, const EcPoint *Q);
-  int (*modOp)(BigNum *, const BigNum *, const BigNum *);
-  int (*checkPubOp)(const struct EcCurve *curve, const EcPoint *pub);
-  int (*randomOp)(const struct EcCurve *curve, BigNum *d);
-  int (*mulScalar)(const struct EcCurve *curve, EcPoint *R, const EcPoint *P, const BigNum *scalar);
+} NistECpCurve;
 
-};
-
-typedef struct EcCurve EcCurve;
-typedef EcCurve NistECpCurve;
+/**
+ * \brief This structure contains the x, y affine coordinates and the z value if we
+ *        use projective coordinates during EC point arithmetic.
+ */
+typedef struct _EcPoint {
+    BigNum *x, *y, *z;
+    BigNum tx, ty, tz;
+} EcPoint;
 
 /**
  * \brief          Marco to initialize a EC point structure.
@@ -112,16 +94,7 @@
  *
  * \param P        Address of the EC point structure.
  */
-#define SET_EC_BASE_POINT(C, P) {EcPoint *e = P;  const EcCurve *c = C; bnCopy(e->x, c->Gx); bnCopy(e->y, c->Gy); bnSetQ(e->z, 1);}
-
-/*
- * EC point helper functions
- */
-extern void ecInitPoint(EcPoint *P);
-
-extern void ecFreePoint(EcPoint *P);
-
-extern void ecSetBasePoint(EcCurve *C, EcPoint *P);
+#define SET_EC_BASE_POINT(C, P) {EcPoint *e = P;  const NistECpCurve *c = C; bnCopy(e->x, c->Gx); bnCopy(e->y, c->Gy); bnSetQ(e->z, 1);}
 
 /**
  * \brief          Get NIST EC curve parameters.
@@ -131,23 +104,23 @@
  *
  * \param curveId  Which curve to initialize
  *
- * \param curve    Pointer to a EcCurve structure
+ * \param curve    Pointer to a NistECpCurve structure
  *
- * \return         0 if successful
+ * \return         0 if successful, or a POLARSSL_ERR_EC_XXX/ POLARSSL_ERR_MPI_XXX error code.
  *
  * \note           Call ecFreeCurveNistECp to return allocated memory.
  */
-int ecGetCurveNistECp(Curves curveId, NistECpCurve *curve);
+int ecGetCurveNistECp(NistCurves curveId, NistECpCurve *curve);
 
 
 /**
- * \brief          Free EC curve parameters.
+ * \brief          Free NIST EC curve parameters.
  *
- * \param curve    Pointer to a EcCurve structure
+ * \param curve    Pointer to a NistECpCurve structure
  *
  * \note           Curve parameters must be initialized calling ecGetCurveNistECp.
  */
-void ecFreeCurveNistECp(EcCurve *curve);
+void ecFreeCurveNistECp(NistECpCurve *curve);
 
 /**
  * \brief          Double an EC point.
@@ -156,13 +129,13 @@
  *                 further reference see RFC 6090 or the standard work <i>Guide to Elliptic
  *                 Curve Cryptography</i>.
  *
- * \param          curve  Address of EC curve structure
+ * \param          curve  Address of Nist EC curve structure
  * \param          R      Address of resulting EC point structure
  * \param          P      Address of the EC point structure
  *
- * \return         0 if successful
+ * \return         0 if successful, or a POLARSSL_ERR_EC_XXX / POLARSSL_ERR_MPI_XXX error code.
  */
-int ecDoublePoint(const EcCurve *curve, EcPoint *R, const EcPoint *P);
+int ecDoublePoint(const NistECpCurve *curve, EcPoint *R, const EcPoint *P);
 
 /**
  * \brief          Add two EC points.
@@ -171,43 +144,43 @@
  *                 further reference see RFC 6090 or the standard work <i>Guide to Elliptic
  *                 Curve Cryptography</i>.
  *
- * \param          curve  Address of EC curve structure
+ * \param          curve  Address of Nist EC curve structure
  * \param          R      Address of resulting EC point structure
  * \param          P      Address of the first EC point structure
  * \param          Q      Address of the second EC point structure
  *
- * \return         0 if successful
+ * \return         0 if successful, or a POLARSSL_ERR_EC_XXX / POLARSSL_ERR_MPI_XXX error code.
  */
-int ecAddPoint(const EcCurve *curve, EcPoint *R, const EcPoint *P, const EcPoint *Q);
+int ecAddPoint(const NistECpCurve *curve, EcPoint *R, const EcPoint *P, const EcPoint *Q);
 
 /**
  * \brief          Mulitply an EC point with a scalar value.
  *
- * \param          curve  Address of EC curve structure
+ * \param          curve  Address of Nist EC curve structure
  * \param          R      Address of resulting EC point structure
  * \param          P      Address of the EC point structure
  * \param          scalar Address of the scalar multi-precision integer value
  *
- * \return         0 if successful
+ * \return         0 if successful, or a POLARSSL_ERR_EC_XXX / POLARSSL_ERR_MPI_XXX error code.
  */
-int ecMulPointScalar(const EcCurve *curve, EcPoint *R, const EcPoint *P, const BigNum *scalar);
+int ecMulPointScalar(const NistECpCurve *curve, EcPoint *R, const EcPoint *P, const BigNum *scalar);
 
 /**
  * \brief          Convert an EC point from Jacobian projective coordinates to normal affine x/y coordinates.
  *
- * \param          curve  Address of EC curve structure
+ * \param          curve  Address of Nist EC curve structure
  * \param          R      Address of EC point structure that receives the x/y coordinates
  * \param          P      Address of the EC point structure that contains the jacobian x/y/z coordinates.
  *
- * \return         0 if successful
+ * \return         0 if successful, or a POLARSSL_ERR_EC_XXX / POLARSSL_ERR_MPI_XXX error code.
  */
-int ecGetAffine(const EcCurve *curve, EcPoint *R, const EcPoint *P);
+int ecGetAffine(const NistECpCurve *curve, EcPoint *R, const EcPoint *P);
 
 /**
  * @brief Generate a random number.
  *
  * The method generates a random number and checks if it matches the curve restricitions.
- * Use this number as ECDH private key.
+ * Use this number to generate a ECDH public key.
  *
  * @param curve the NIST curve to use.
  *
@@ -215,35 +188,6 @@
  */
 int ecGenerateRandomNumber(const NistECpCurve *curve, BigNum *d);
 
-/**
- * @brief Check a public key.
- *
- * The method checks if a public key is valid. For NIST curves it uses the
- * ECC Partial Validation, NIST SP800-56A section 5.6.2.6
- * 
- * For other curves it computes the equation and compares the left hand and 
- * the right handresults. If they are equal the point is on the curve.
- *
- * @param curve the curve to use.
- *
- * @param pub the public key to check.
- *
- * @returns true (!0) if the check was ok, false (0) otherwise.
- *
- * @note The function uses some scratch pad variable of the NistECpCurve structure.
- */
-int ecCheckPubKey(const EcCurve *curve, const EcPoint *pub);
-
-int ecGetCurvesCurve(Curves curveId, EcCurve *curve);
-
-void ecFreeCurvesCurve(EcCurve *curve);
-
-/**
- * This is a special function for DJB's curve 25519. Actually it's the scalar multiplication
- * mypublic = basepoint * secret
- */
-int curve25519_donna(unsigned char *mypublic, const unsigned char *secret, const unsigned char *basepoint);
-
 /*
  * Some additional functions that are not available in bnlib
  */
@@ -255,11 +199,11 @@
 
 int bnSubQMod_ (struct BigNum *rslt, unsigned n1, struct BigNum *mod);
 
-int bnMulMod_ (struct BigNum *rslt, struct BigNum *n1, struct BigNum *n2, struct BigNum *mod, const EcCurve *curve);
+int bnMulMod_ (struct BigNum *rslt, struct BigNum *n1, struct BigNum *n2, struct BigNum *mod);
 
-int bnMulQMod_ (struct BigNum *rslt, struct BigNum *n1, unsigned n2, struct BigNum *mod, const EcCurve *curve);
+int bnMulQMod_ (struct BigNum *rslt, struct BigNum *n1, unsigned n2, struct BigNum *mod);
 
-int bnSquareMod_ (struct BigNum *rslt, struct BigNum *n1, struct BigNum *mod, const EcCurve *curve);
+int bnSquareMod_ (struct BigNum *rslt, struct BigNum *n1, struct BigNum *mod);
 
 #ifdef __cplusplus
 }
diff --git a/jni/libzrtp/sources/bnlib/ec/ecdh.c b/jni/libzrtp/sources/bnlib/ec/ecdh.c
index 8d1bc23..25ea8cd 100644
--- a/jni/libzrtp/sources/bnlib/ec/ecdh.c
+++ b/jni/libzrtp/sources/bnlib/ec/ecdh.c
@@ -9,7 +9,7 @@
 #include <ec/ec.h>
 #include <ec/ecdh.h>
 
-int ecdhGeneratePublic(const EcCurve *curve, EcPoint *Q, const BigNum *d)
+int ecdhGeneratePublic(const NistECpCurve *curve, EcPoint *Q, const BigNum *d)
 {
     EcPoint G;
 
@@ -21,10 +21,10 @@
 
     FREE_EC_POINT(&G);
 
-    return ecCheckPubKey(curve, Q);
+    return 0;
 }
 
-int ecdhComputeAgreement(const EcCurve *curve, BigNum *agreement, const EcPoint *Q, const BigNum *d)
+int ecdhComputeAgreement(const NistECpCurve *curve, BigNum *agreement, const EcPoint *Q, const BigNum *d)
 {
     EcPoint t0;
 
diff --git a/jni/libzrtp/sources/bnlib/ec/ecdh.h b/jni/libzrtp/sources/bnlib/ec/ecdh.h
index 7ec32ad..0cf083a 100644
--- a/jni/libzrtp/sources/bnlib/ec/ecdh.h
+++ b/jni/libzrtp/sources/bnlib/ec/ecdh.h
@@ -22,17 +22,15 @@
 /**
  * @brief Takes a secret large random number and computes the public EC point.
  *
- * @param curve is the curve to use.
+ * @param curve is the NIST curve to use.
  *
  * @param Q the functions writes the computed public point in this parameter.
  *
  * @param d is the secret random number.
  *
- * @return @c true (!0) if public key was computed, @c false otherwise.
- *
  * @sa ecGenerateRandomNumber
  */
-int ecdhGeneratePublic(const EcCurve *curve, EcPoint *Q, const BigNum *d);
+int ecdhGeneratePublic(const NistECpCurve *curve, EcPoint *Q, const BigNum *d);
 
 /**
  * @brief Computes the key agreement value.
@@ -40,7 +38,7 @@
  * Takes the public EC point of the other party and applies the EC DH algorithm
  * to compute the agreed value.
  *
- * @param curve is the curve to use, must be the same curve as used in
+ * @param curve is the NIST curve to use, must be the same curve as used in
  *              @c ecdhGeneratePublic.
  *
  * @param agreemtn the functions writes the computed agreed value in this parameter.
@@ -49,7 +47,7 @@
  *
  * @param d is the secret random number.
  */
-int ecdhComputeAgreement(const EcCurve *curve, BigNum *agreement, const EcPoint *Q, const BigNum *d);
+int ecdhComputeAgreement(const NistECpCurve *curve, BigNum *agreement, const EcPoint *Q, const BigNum *d);
 
 #ifdef __cplusplus
 }
diff --git a/jni/libzrtp/sources/bnlib/lbn32.c b/jni/libzrtp/sources/bnlib/lbn32.c
index 73fedcb..831b25c 100644
--- a/jni/libzrtp/sources/bnlib/lbn32.c
+++ b/jni/libzrtp/sources/bnlib/lbn32.c
@@ -69,7 +69,7 @@
 #define HAVE_CONFIG_H 0
 #endif
 #if HAVE_CONFIG_H
-#include <bnconfig.h>
+#include "bnconfig.h"
 #endif
 
 /*
diff --git a/jni/libzrtp/sources/bnlib/lbn64.c b/jni/libzrtp/sources/bnlib/lbn64.c
index e930652..f14deea 100644
--- a/jni/libzrtp/sources/bnlib/lbn64.c
+++ b/jni/libzrtp/sources/bnlib/lbn64.c
@@ -69,7 +69,7 @@
 #define HAVE_CONFIG_H 0
 #endif
 #if HAVE_CONFIG_H
-#include <bnconfig.h>
+#include "bnconfig.h"
 #endif
 
 /*