diff --git a/.gitignore b/.gitignore index 56d45d7..5d60f4f 100644 --- a/.gitignore +++ b/.gitignore @@ -7,6 +7,7 @@ *.prof *.out *.log +src/*_tune.c tmp/ doc/ inc/blst.h* diff --git a/src/Makefile b/src/Makefile index 3817959..5fb5daf 100644 --- a/src/Makefile +++ b/src/Makefile @@ -16,14 +16,14 @@ INCLUDE_DIRS = ../inc libckzg.a: $(LIB_OBJ) Makefile ar rc libckzg.a $(LIB_OBJ) -%_test: %_test.c debug_util.o test_util.o libckzg.a Makefile - clang -Wall -I$(INCLUDE_DIRS) $(CFLAGS) $(KZG_CFLAGS) -o $@ $@.c debug_util.o test_util.o libckzg.a -L../lib -lblst +%_test: %.c debug_util.o test_util.o libckzg.a Makefile + clang -Wall -DKZGTEST -I$(INCLUDE_DIRS) $(CFLAGS) $(KZG_CFLAGS) -o $@ $*.c debug_util.o test_util.o libckzg.a -L../lib -lblst ./$@ # This version will abort on error and print the file and line number -%_test_debug: KZG_CFLAGS += -g -O0 -DDEBUG -%_test_debug: %_test.c debug_util.o test_util.o libckzg.a Makefile - clang -Wall -I$(INCLUDE_DIRS) $(CFLAGS) $(KZG_CFLAGS) -o $@ $*_test.c debug_util.o test_util.o libckzg.a -L../lib -lblst +%_test_debug: KZG_CFLAGS += -g -O0 +%_test_debug: %.c debug_util.o test_util.o libckzg.a Makefile + clang -Wall -DKZGTEST -DDEBUG -I$(INCLUDE_DIRS) $(CFLAGS) $(KZG_CFLAGS) -o $@ $*.c debug_util.o test_util.o libckzg.a -L../lib -lblst # Benchmarks %_bench: KZG_CFLAGS += -O diff --git a/src/bls12_381.c b/src/bls12_381.c index 3f48a0d..2b758cf 100644 --- a/src/bls12_381.c +++ b/src/bls12_381.c @@ -35,7 +35,7 @@ * @param[in] b A non-zero byte * @return The index of the highest set bit */ -int log_2_byte(byte b) { +static int log_2_byte(byte b) { int r, shift; r = (b > 0xF) << 2; b >>= r; @@ -472,4 +472,256 @@ bool pairings_verify(const g1_t *a1, const g2_t *a2, const g1_t *b1, const g2_t return blst_fp12_is_one(>_point); } -#endif // BLST \ No newline at end of file +#endif // BLST + +#ifdef KZGTEST + +#include "../inc/acutest.h" +#include "test_util.h" + +// This is -1 (the second root of unity) +uint64_t m1[] = {0xffffffff00000000L, 0x53bda402fffe5bfeL, 0x3339d80809a1d805L, 0x73eda753299d7d48L}; + +void log_2_byte_works(void) { + // TEST_CHECK(0 == log_2_byte(0x00)); + TEST_CHECK(0 == log_2_byte(0x01)); + TEST_CHECK(7 == log_2_byte(0x80)); + TEST_CHECK(7 == log_2_byte(0xff)); + TEST_CHECK(4 == log_2_byte(0x10)); +} + +void fr_is_zero_works(void) { + fr_t zero; + fr_from_uint64(&zero, 0); + TEST_CHECK(fr_is_zero(&zero)); +} + +void fr_is_one_works(void) { + TEST_CHECK(fr_is_one(&fr_one)); +} + +void fr_is_null_works(void) { + TEST_CHECK(fr_is_null(&fr_null)); + TEST_CHECK(!fr_is_null(&fr_zero)); + TEST_CHECK(!fr_is_null(&fr_one)); +} + +void fr_from_uint64_works(void) { + fr_t a; + fr_from_uint64(&a, 1); + TEST_CHECK(fr_is_one(&a)); +} + +void fr_equal_works(void) { + // A couple of arbitrary roots of unity + uint64_t aa[] = {0x0001000000000000L, 0xec03000276030000L, 0x8d51ccce760304d0L, 0x0000000000000000L}; + uint64_t bb[] = {0x8dd702cb688bc087L, 0xa032824078eaa4feL, 0xa733b23a98ca5b22L, 0x3f96405d25a31660L}; + fr_t a, b; + fr_from_uint64s(&a, aa); + fr_from_uint64s(&b, bb); + TEST_CHECK(true == fr_equal(&a, &a)); + TEST_CHECK(false == fr_equal(&a, &b)); +} + +void fr_negate_works(void) { + fr_t minus1, res; + fr_from_uint64s(&minus1, m1); + fr_negate(&res, &minus1); + TEST_CHECK(fr_is_one(&res)); +} + +void fr_pow_works(void) { + // a^pow + uint64_t pow = 123456; + fr_t a, expected, actual; + fr_from_uint64(&a, 197); + + // Do it the slow way + expected = fr_one; + for (uint64_t i = 0; i < pow; i++) { + fr_mul(&expected, &expected, &a); + } + + // Do it the quick way + fr_pow(&actual, &a, pow); + + TEST_CHECK(fr_equal(&expected, &actual)); +} + +void fr_div_works(void) { + fr_t a, b, tmp, actual; + + fr_from_uint64(&a, 197); + fr_from_uint64(&b, 123456); + + fr_div(&tmp, &a, &b); + fr_mul(&actual, &tmp, &b); + + TEST_CHECK(fr_equal(&a, &actual)); +} + +// This is strictly undefined, but conventionally 0 is returned +void fr_div_by_zero(void) { + fr_t a, b, tmp; + + fr_from_uint64(&a, 197); + fr_from_uint64(&b, 0); + + fr_div(&tmp, &a, &b); + + TEST_CHECK(fr_is_zero(&tmp)); +} + +void fr_uint64s_roundtrip(void) { + fr_t fr; + uint64_t expected[4] = {1, 2, 3, 4}; + uint64_t actual[4]; + + fr_from_uint64s(&fr, expected); + fr_to_uint64s(actual, &fr); + + TEST_CHECK(expected[0] == actual[0]); + TEST_CHECK(expected[1] == actual[1]); + TEST_CHECK(expected[2] == actual[2]); + TEST_CHECK(expected[3] == actual[3]); +} + +void p1_mul_works(void) { + fr_t minus1; + g1_t res; + + // Multiply the generator by minus one (the second root of unity) + fr_from_uint64s(&minus1, m1); + g1_mul(&res, &g1_generator, &minus1); + + // We should end up with negative the generator + TEST_CHECK(g1_equal(&res, &g1_negative_generator)); +} + +void p1_sub_works(void) { + g1_t tmp, res; + + // 2 * g1_gen = g1_gen - g1_gen_neg + g1_dbl(&tmp, &g1_generator); + g1_sub(&res, &g1_generator, &g1_negative_generator); + + TEST_CHECK(g1_equal(&tmp, &res)); +} + +void p2_mul_works(void) { + fr_t minus1; + g2_t res; + + // Multiply the generator by minus one (the second root of unity) + fr_from_uint64s(&minus1, m1); + g2_mul(&res, &g2_generator, &minus1); + + TEST_CHECK(g2_equal(&res, &g2_negative_generator)); +} + +void p2_sub_works(void) { + g2_t tmp, res; + + // 2 * g2_gen = g2_gen - g2_gen_neg + g2_dbl(&tmp, &g2_generator); + g2_sub(&res, &g2_generator, &g2_negative_generator); + + TEST_CHECK(g2_equal(&tmp, &res)); +} + +void g1_identity_is_infinity(void) { + TEST_CHECK(g1_is_inf(&g1_identity)); +} + +void g1_identity_is_identity(void) { + g1_t actual; + g1_add_or_dbl(&actual, &g1_generator, &g1_identity); + TEST_CHECK(g1_equal(&g1_generator, &actual)); +} + +void g1_make_linear_combination(void) { + int len = 255; + fr_t coeffs[len], tmp; + g1_t p[len], res, exp; + for (int i = 0; i < len; i++) { + fr_from_uint64(coeffs + i, i + 1); + p[i] = g1_generator; + } + + // Expected result + fr_from_uint64(&tmp, len * (len + 1) / 2); + g1_mul(&exp, &g1_generator, &tmp); + + // Test result + g1_linear_combination(&res, p, coeffs, len); + TEST_CHECK(g1_equal(&exp, &res)); +} + +void g1_random_linear_combination(void) { + int len = 8192; + fr_t coeffs[len]; + g1_t p[len], p1tmp = g1_generator; + for (int i = 0; i < len; i++) { + coeffs[i] = rand_fr(); + p[i] = p1tmp; + g1_dbl(&p1tmp, &p1tmp); + } + + // Expected result + g1_t exp = g1_identity; + for (uint64_t i = 0; i < len; i++) { + g1_mul(&p1tmp, &p[i], &coeffs[i]); + g1_add_or_dbl(&exp, &exp, &p1tmp); + } + + // Test result + g1_t res; + g1_linear_combination(&res, p, coeffs, len); + TEST_CHECK(g1_equal(&exp, &res)); +} + +void pairings_work(void) { + // Verify that e([3]g1, [5]g2) = e([5]g1, [3]g2) + fr_t three, five; + g1_t g1_3, g1_5; + g2_t g2_3, g2_5; + + // Set up + fr_from_uint64(&three, 3); + fr_from_uint64(&five, 5); + g1_mul(&g1_3, &g1_generator, &three); + g1_mul(&g1_5, &g1_generator, &five); + g2_mul(&g2_3, &g2_generator, &three); + g2_mul(&g2_5, &g2_generator, &five); + + // Verify the pairing + TEST_CHECK(true == pairings_verify(&g1_3, &g2_5, &g1_5, &g2_3)); + TEST_CHECK(false == pairings_verify(&g1_3, &g2_3, &g1_5, &g2_5)); +} + +TEST_LIST = { + {"BLS12_384_TEST", title}, + {"log_2_byte_works", log_2_byte_works}, + {"fr_is_zero_works", fr_is_zero_works}, + {"fr_is_one_works", fr_is_one_works}, + {"fr_is_null_works", fr_is_null_works}, + {"fr_from_uint64_works", fr_from_uint64_works}, + {"fr_equal_works", fr_equal_works}, + {"fr_negate_works", fr_negate_works}, + {"fr_pow_works", fr_pow_works}, + {"fr_div_works", fr_div_works}, + {"fr_div_by_zero", fr_div_by_zero}, + {"fr_uint64s_roundtrip", fr_uint64s_roundtrip}, + {"p1_mul_works", p1_mul_works}, + {"p1_sub_works", p1_sub_works}, + {"p2_mul_works", p2_mul_works}, + {"p2_sub_works", p2_sub_works}, + {"g1_identity_is_infinity", g1_identity_is_infinity}, + {"g1_identity_is_identity", g1_identity_is_identity}, + {"g1_make_linear_combination", g1_make_linear_combination}, + {"g1_random_linear_combination", g1_make_linear_combination}, + {"pairings_work", pairings_work}, + {NULL, NULL} /* zero record marks the end of the list */ +}; + +#endif // KZGTEST \ No newline at end of file diff --git a/src/bls12_381.h b/src/bls12_381.h index ed6e44c..7e8a236 100644 --- a/src/bls12_381.h +++ b/src/bls12_381.h @@ -97,7 +97,6 @@ static const g2_t g2_negative_generator = {{{{0xf5f28fa202940a10L, 0xb3f5fb2687b #endif // BLST // All the functions in the interface -int log_2_byte(byte b); bool fr_is_zero(const fr_t *p); bool fr_is_one(const fr_t *p); bool fr_is_null(const fr_t *p); diff --git a/src/bls12_381_test.c b/src/bls12_381_test.c deleted file mode 100644 index 35216cf..0000000 --- a/src/bls12_381_test.c +++ /dev/null @@ -1,263 +0,0 @@ -/* - * Copyright 2021 Benjamin Edgington - * - * Licensed under the Apache License, Version 2.0 (the "License"); - * you may not use this file except in compliance with the License. - * You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -#include "../inc/acutest.h" -#include "test_util.h" - -// This is -1 (the second root of unity) -uint64_t m1[] = {0xffffffff00000000L, 0x53bda402fffe5bfeL, 0x3339d80809a1d805L, 0x73eda753299d7d48L}; - -void log_2_byte_works(void) { - // TEST_CHECK(0 == log_2_byte(0x00)); - TEST_CHECK(0 == log_2_byte(0x01)); - TEST_CHECK(7 == log_2_byte(0x80)); - TEST_CHECK(7 == log_2_byte(0xff)); - TEST_CHECK(4 == log_2_byte(0x10)); -} - -void fr_is_zero_works(void) { - fr_t zero; - fr_from_uint64(&zero, 0); - TEST_CHECK(fr_is_zero(&zero)); -} - -void fr_is_one_works(void) { - TEST_CHECK(fr_is_one(&fr_one)); -} - -void fr_is_null_works(void) { - TEST_CHECK(fr_is_null(&fr_null)); - TEST_CHECK(!fr_is_null(&fr_zero)); - TEST_CHECK(!fr_is_null(&fr_one)); -} - -void fr_from_uint64_works(void) { - fr_t a; - fr_from_uint64(&a, 1); - TEST_CHECK(fr_is_one(&a)); -} - -void fr_equal_works(void) { - // A couple of arbitrary roots of unity - uint64_t aa[] = {0x0001000000000000L, 0xec03000276030000L, 0x8d51ccce760304d0L, 0x0000000000000000L}; - uint64_t bb[] = {0x8dd702cb688bc087L, 0xa032824078eaa4feL, 0xa733b23a98ca5b22L, 0x3f96405d25a31660L}; - fr_t a, b; - fr_from_uint64s(&a, aa); - fr_from_uint64s(&b, bb); - TEST_CHECK(true == fr_equal(&a, &a)); - TEST_CHECK(false == fr_equal(&a, &b)); -} - -void fr_negate_works(void) { - fr_t minus1, res; - fr_from_uint64s(&minus1, m1); - fr_negate(&res, &minus1); - TEST_CHECK(fr_is_one(&res)); -} - -void fr_pow_works(void) { - // a^pow - uint64_t pow = 123456; - fr_t a, expected, actual; - fr_from_uint64(&a, 197); - - // Do it the slow way - expected = fr_one; - for (uint64_t i = 0; i < pow; i++) { - fr_mul(&expected, &expected, &a); - } - - // Do it the quick way - fr_pow(&actual, &a, pow); - - TEST_CHECK(fr_equal(&expected, &actual)); -} - -void fr_div_works(void) { - fr_t a, b, tmp, actual; - - fr_from_uint64(&a, 197); - fr_from_uint64(&b, 123456); - - fr_div(&tmp, &a, &b); - fr_mul(&actual, &tmp, &b); - - TEST_CHECK(fr_equal(&a, &actual)); -} - -// This is strictly undefined, but conventionally 0 is returned -void fr_div_by_zero(void) { - fr_t a, b, tmp; - - fr_from_uint64(&a, 197); - fr_from_uint64(&b, 0); - - fr_div(&tmp, &a, &b); - - TEST_CHECK(fr_is_zero(&tmp)); -} - -void fr_uint64s_roundtrip(void) { - fr_t fr; - uint64_t expected[4] = {1, 2, 3, 4}; - uint64_t actual[4]; - - fr_from_uint64s(&fr, expected); - fr_to_uint64s(actual, &fr); - - TEST_CHECK(expected[0] == actual[0]); - TEST_CHECK(expected[1] == actual[1]); - TEST_CHECK(expected[2] == actual[2]); - TEST_CHECK(expected[3] == actual[3]); -} - -void p1_mul_works(void) { - fr_t minus1; - g1_t res; - - // Multiply the generator by minus one (the second root of unity) - fr_from_uint64s(&minus1, m1); - g1_mul(&res, &g1_generator, &minus1); - - // We should end up with negative the generator - TEST_CHECK(g1_equal(&res, &g1_negative_generator)); -} - -void p1_sub_works(void) { - g1_t tmp, res; - - // 2 * g1_gen = g1_gen - g1_gen_neg - g1_dbl(&tmp, &g1_generator); - g1_sub(&res, &g1_generator, &g1_negative_generator); - - TEST_CHECK(g1_equal(&tmp, &res)); -} - -void p2_mul_works(void) { - fr_t minus1; - g2_t res; - - // Multiply the generator by minus one (the second root of unity) - fr_from_uint64s(&minus1, m1); - g2_mul(&res, &g2_generator, &minus1); - - TEST_CHECK(g2_equal(&res, &g2_negative_generator)); -} - -void p2_sub_works(void) { - g2_t tmp, res; - - // 2 * g2_gen = g2_gen - g2_gen_neg - g2_dbl(&tmp, &g2_generator); - g2_sub(&res, &g2_generator, &g2_negative_generator); - - TEST_CHECK(g2_equal(&tmp, &res)); -} - -void g1_identity_is_infinity(void) { - TEST_CHECK(g1_is_inf(&g1_identity)); -} - -void g1_identity_is_identity(void) { - g1_t actual; - g1_add_or_dbl(&actual, &g1_generator, &g1_identity); - TEST_CHECK(g1_equal(&g1_generator, &actual)); -} - -void g1_make_linear_combination(void) { - int len = 255; - fr_t coeffs[len], tmp; - g1_t p[len], res, exp; - for (int i = 0; i < len; i++) { - fr_from_uint64(coeffs + i, i + 1); - p[i] = g1_generator; - } - - // Expected result - fr_from_uint64(&tmp, len * (len + 1) / 2); - g1_mul(&exp, &g1_generator, &tmp); - - // Test result - g1_linear_combination(&res, p, coeffs, len); - TEST_CHECK(g1_equal(&exp, &res)); -} - -void g1_random_linear_combination(void) { - int len = 8192; - fr_t coeffs[len]; - g1_t p[len], p1tmp = g1_generator; - for (int i = 0; i < len; i++) { - coeffs[i] = rand_fr(); - p[i] = p1tmp; - blst_p1_double(&p1tmp, &p1tmp); - } - - // Expected result - g1_t exp = g1_identity; - for (uint64_t i = 0; i < len; i++) { - g1_mul(&p1tmp, &p[i], &coeffs[i]); - blst_p1_add_or_double(&exp, &exp, &p1tmp); - } - - // Test result - g1_t res; - g1_linear_combination(&res, p, coeffs, len); - TEST_CHECK(g1_equal(&exp, &res)); -} - -void pairings_work(void) { - // Verify that e([3]g1, [5]g2) = e([5]g1, [3]g2) - fr_t three, five; - g1_t g1_3, g1_5; - g2_t g2_3, g2_5; - - // Set up - fr_from_uint64(&three, 3); - fr_from_uint64(&five, 5); - g1_mul(&g1_3, &g1_generator, &three); - g1_mul(&g1_5, &g1_generator, &five); - g2_mul(&g2_3, &g2_generator, &three); - g2_mul(&g2_5, &g2_generator, &five); - - // Verify the pairing - TEST_CHECK(true == pairings_verify(&g1_3, &g2_5, &g1_5, &g2_3)); - TEST_CHECK(false == pairings_verify(&g1_3, &g2_3, &g1_5, &g2_5)); -} - -TEST_LIST = { - {"BLS12_384_TEST", title}, - {"log_2_byte_works", log_2_byte_works}, - {"fr_is_zero_works", fr_is_zero_works}, - {"fr_is_one_works", fr_is_one_works}, - {"fr_is_null_works", fr_is_null_works}, - {"fr_from_uint64_works", fr_from_uint64_works}, - {"fr_equal_works", fr_equal_works}, - {"fr_negate_works", fr_negate_works}, - {"fr_pow_works", fr_pow_works}, - {"fr_div_works", fr_div_works}, - {"fr_div_by_zero", fr_div_by_zero}, - {"fr_uint64s_roundtrip", fr_uint64s_roundtrip}, - {"p1_mul_works", p1_mul_works}, - {"p1_sub_works", p1_sub_works}, - {"p2_mul_works", p2_mul_works}, - {"p2_sub_works", p2_sub_works}, - {"g1_identity_is_infinity", g1_identity_is_infinity}, - {"g1_identity_is_identity", g1_identity_is_identity}, - {"g1_make_linear_combination", g1_make_linear_combination}, - {"g1_random_linear_combination", g1_make_linear_combination}, - {"pairings_work", pairings_work}, - {NULL, NULL} /* zero record marks the end of the list */ -}; diff --git a/src/c_kzg_util.c b/src/c_kzg_util.c index cfd5a93..31fcf81 100644 --- a/src/c_kzg_util.c +++ b/src/c_kzg_util.c @@ -30,7 +30,7 @@ * @retval C_CZK_OK All is well * @retval C_CZK_MALLOC Memory allocation failed */ -C_KZG_RET c_kzg_malloc(void **x, size_t n) { +static C_KZG_RET c_kzg_malloc(void **x, size_t n) { if (n > 0) { *x = malloc(n); return *x != NULL ? C_KZG_OK : C_KZG_MALLOC; @@ -136,3 +136,28 @@ C_KZG_RET new_g2_array(g2_t **x, size_t n) { C_KZG_RET new_poly_array(poly **x, size_t n) { return c_kzg_malloc((void **)x, n * sizeof **x); } + +#ifdef KZGTEST + +#include "../inc/acutest.h" +#include "test_util.h" + +void malloc_works(void) { + int *p; + TEST_CHECK(C_KZG_OK == c_kzg_malloc((void **)&p, 4)); + free(p); +} + +void malloc_huge_fails(void) { + int *p; + TEST_CHECK(C_KZG_MALLOC == c_kzg_malloc((void **)&p, -1)); +} + +TEST_LIST = { + {"C_KZG_UTIL_TEST", title}, + {"malloc_works", malloc_works}, + {"malloc_huge_fails", malloc_huge_fails}, + {NULL, NULL} /* zero record marks the end of the list */ +}; + +#endif // KZGTEST \ No newline at end of file diff --git a/src/c_kzg_util.h b/src/c_kzg_util.h index d6230cc..72603d4 100644 --- a/src/c_kzg_util.h +++ b/src/c_kzg_util.h @@ -42,7 +42,6 @@ _a > _b ? _a : _b; \ }) -C_KZG_RET c_kzg_malloc(void **p, size_t n); C_KZG_RET new_uint64_array(uint64_t **x, size_t n); C_KZG_RET new_fr_array(fr_t **x, size_t n); C_KZG_RET new_fr_array_2(fr_t ***x, size_t n); diff --git a/src/c_kzg_util_test.c b/src/c_kzg_util_test.c deleted file mode 100644 index 089f676..0000000 --- a/src/c_kzg_util_test.c +++ /dev/null @@ -1,37 +0,0 @@ -/* - * Copyright 2021 Benjamin Edgington - * - * Licensed under the Apache License, Version 2.0 (the "License"); - * you may not use this file except in compliance with the License. - * You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -#include "../inc/acutest.h" -#include "test_util.h" -#include "c_kzg_util.h" - -void malloc_works(void) { - int *p; - TEST_CHECK(C_KZG_OK == c_kzg_malloc((void **)&p, 4)); - free(p); -} - -void malloc_huge_fails(void) { - int *p; - TEST_CHECK(C_KZG_MALLOC == c_kzg_malloc((void **)&p, -1)); -} - -TEST_LIST = { - {"C_KZG_UTIL_TEST", title}, - {"malloc_works", malloc_works}, - {"malloc_huge_fails", malloc_huge_fails}, - {NULL, NULL} /* zero record marks the end of the list */ -}; diff --git a/src/das_extension.c b/src/das_extension.c index 92b3f91..ce1eda1 100644 --- a/src/das_extension.c +++ b/src/das_extension.c @@ -110,4 +110,101 @@ C_KZG_RET das_fft_extension(fr_t *vals, uint64_t n, const FFTSettings *fs) { } return C_KZG_OK; -} \ No newline at end of file +} + +#ifdef KZGTEST + +#include "../inc/acutest.h" +#include "c_kzg_util.h" +#include "test_util.h" + +void das_extension_test_known(void) { + FFTSettings fs; + uint64_t half; + fr_t *data; + // The expected output (from go-kzg): + const uint64_t expected_u[8][4] = { + {0xa0c43757db972d7dL, 0x79d15a1e0677962cL, 0xf678865c0c95fa6aL, 0x4e85fd4814f96825L}, + {0xad9f844939f2705dL, 0x319e440c9f3b0325L, 0x4cbd29a60e160a28L, 0x665961d85d90c4c0L}, + {0x5f3ac8a72468d28bL, 0xede949e28383c5d2L, 0xaf6f84dd8708d8c9L, 0x2567aa0b14a41521L}, + {0x25abe312b96aadadL, 0x4abf043f091ff417L, 0x43824b53e09536dbL, 0x195dbe06a28ca227L}, + {0x5f3ac8a72468d28bL, 0xede949e28383c5d2L, 0xaf6f84dd8708d8c9L, 0x2567aa0b14a41521L}, + {0xad9f844939f2705dL, 0x319e440c9f3b0325L, 0x4cbd29a60e160a28L, 0x665961d85d90c4c0L}, + {0xa0c43757db972d7dL, 0x79d15a1e0677962cL, 0xf678865c0c95fa6aL, 0x4e85fd4814f96825L}, + {0x7f171458d2b071a9L, 0xd185bbb2a46cbd9bL, 0xa41aab0d02886e80L, 0x01cacceef58ccee9L}}; + + TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 4)); + half = fs.max_width / 2; + + TEST_CHECK(C_KZG_OK == new_fr_array(&data, half)); + for (uint64_t i = 0; i < half; i++) { + fr_from_uint64(data + i, i); + } + + TEST_CHECK(C_KZG_OK == das_fft_extension(data, half, &fs)); + + // Check against the expected values + for (uint64_t i = 0; i < 8; i++) { + fr_t expected; + fr_from_uint64s(&expected, expected_u[i]); + TEST_CHECK(fr_equal(&expected, data + i)); + } + + free(data); + free_fft_settings(&fs); +} + +void das_extension_test_random(void) { + FFTSettings fs; + fr_t *even_data, *odd_data, *data, *coeffs; + int max_scale = 15; + + TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, max_scale)); + for (int scale = 1; scale <= max_scale; scale++) { + uint64_t width = (uint64_t)1 << scale; + TEST_ASSERT(width <= fs.max_width); + TEST_CHECK(C_KZG_OK == new_fr_array(&even_data, width / 2)); + TEST_CHECK(C_KZG_OK == new_fr_array(&odd_data, width / 2)); + TEST_CHECK(C_KZG_OK == new_fr_array(&data, width)); + TEST_CHECK(C_KZG_OK == new_fr_array(&coeffs, width)); + + for (int rep = 0; rep < 4; rep++) { + + // Make random even data and duplicate temporarily in the odd_data + for (int i = 0; i < width / 2; i++) { + even_data[i] = rand_fr(); + odd_data[i] = even_data[i]; + } + + // Extend the even data to create the odd data required to make the second half of the FFT zero + TEST_CHECK(C_KZG_OK == das_fft_extension(odd_data, width / 2, &fs)); + + // Reconstruct the full data + for (int i = 0; i < width; i += 2) { + data[i] = even_data[i / 2]; + data[i + 1] = odd_data[i / 2]; + } + TEST_CHECK(C_KZG_OK == fft_fr(coeffs, data, true, width, &fs)); + + // Second half of the coefficients should be all zeros + for (int i = width / 2; i < width; i++) { + TEST_CHECK(fr_is_zero(&coeffs[i])); + } + } + + free(even_data); + free(odd_data); + free(data); + free(coeffs); + } + free_fft_settings(&fs); +} + +TEST_LIST = { + {"DAS_EXTENSION_TEST", title}, + {"das_extension_test_known", das_extension_test_known}, + {"das_extension_test_random", das_extension_test_random}, + {NULL, NULL} /* zero record marks the end of the list */ +}; + +#endif // KZGTEST \ No newline at end of file diff --git a/src/das_extension_test.c b/src/das_extension_test.c deleted file mode 100644 index 07e052e..0000000 --- a/src/das_extension_test.c +++ /dev/null @@ -1,110 +0,0 @@ -/* - * Copyright 2021 Benjamin Edgington - * - * Licensed under the Apache License, Version 2.0 (the "License"); - * you may not use this file except in compliance with the License. - * You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -#include "../inc/acutest.h" -#include "c_kzg_util.h" -#include "test_util.h" -#include "fft_fr.h" -#include "das_extension.h" - -void das_extension_test_known(void) { - FFTSettings fs; - uint64_t half; - fr_t *data; - // The expected output (from go-kzg): - const uint64_t expected_u[8][4] = { - {0xa0c43757db972d7dL, 0x79d15a1e0677962cL, 0xf678865c0c95fa6aL, 0x4e85fd4814f96825L}, - {0xad9f844939f2705dL, 0x319e440c9f3b0325L, 0x4cbd29a60e160a28L, 0x665961d85d90c4c0L}, - {0x5f3ac8a72468d28bL, 0xede949e28383c5d2L, 0xaf6f84dd8708d8c9L, 0x2567aa0b14a41521L}, - {0x25abe312b96aadadL, 0x4abf043f091ff417L, 0x43824b53e09536dbL, 0x195dbe06a28ca227L}, - {0x5f3ac8a72468d28bL, 0xede949e28383c5d2L, 0xaf6f84dd8708d8c9L, 0x2567aa0b14a41521L}, - {0xad9f844939f2705dL, 0x319e440c9f3b0325L, 0x4cbd29a60e160a28L, 0x665961d85d90c4c0L}, - {0xa0c43757db972d7dL, 0x79d15a1e0677962cL, 0xf678865c0c95fa6aL, 0x4e85fd4814f96825L}, - {0x7f171458d2b071a9L, 0xd185bbb2a46cbd9bL, 0xa41aab0d02886e80L, 0x01cacceef58ccee9L}}; - - TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 4)); - half = fs.max_width / 2; - - TEST_CHECK(C_KZG_OK == new_fr_array(&data, half)); - for (uint64_t i = 0; i < half; i++) { - fr_from_uint64(data + i, i); - } - - TEST_CHECK(C_KZG_OK == das_fft_extension(data, half, &fs)); - - // Check against the expected values - for (uint64_t i = 0; i < 8; i++) { - fr_t expected; - fr_from_uint64s(&expected, expected_u[i]); - TEST_CHECK(fr_equal(&expected, data + i)); - } - - free(data); - free_fft_settings(&fs); -} - -void das_extension_test_random(void) { - FFTSettings fs; - fr_t *even_data, *odd_data, *data, *coeffs; - int max_scale = 15; - - TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, max_scale)); - for (int scale = 1; scale <= max_scale; scale++) { - uint64_t width = (uint64_t)1 << scale; - TEST_ASSERT(width <= fs.max_width); - TEST_CHECK(C_KZG_OK == new_fr_array(&even_data, width / 2)); - TEST_CHECK(C_KZG_OK == new_fr_array(&odd_data, width / 2)); - TEST_CHECK(C_KZG_OK == new_fr_array(&data, width)); - TEST_CHECK(C_KZG_OK == new_fr_array(&coeffs, width)); - - for (int rep = 0; rep < 4; rep++) { - - // Make random even data and duplicate temporarily in the odd_data - for (int i = 0; i < width / 2; i++) { - even_data[i] = rand_fr(); - odd_data[i] = even_data[i]; - } - - // Extend the even data to create the odd data required to make the second half of the FFT zero - TEST_CHECK(C_KZG_OK == das_fft_extension(odd_data, width / 2, &fs)); - - // Reconstruct the full data - for (int i = 0; i < width; i += 2) { - data[i] = even_data[i / 2]; - data[i + 1] = odd_data[i / 2]; - } - TEST_CHECK(C_KZG_OK == fft_fr(coeffs, data, true, width, &fs)); - - // Second half of the coefficients should be all zeros - for (int i = width / 2; i < width; i++) { - TEST_CHECK(fr_is_zero(&coeffs[i])); - } - } - - free(even_data); - free(odd_data); - free(data); - free(coeffs); - } - free_fft_settings(&fs); -} - -TEST_LIST = { - {"DAS_EXTENSION_TEST", title}, - {"das_extension_test_known", das_extension_test_known}, - {"das_extension_test_random", das_extension_test_random}, - {NULL, NULL} /* zero record marks the end of the list */ -}; diff --git a/src/fft_common.c b/src/fft_common.c index e0f0fee..7ec8342 100644 --- a/src/fft_common.c +++ b/src/fft_common.c @@ -35,7 +35,7 @@ * @retval C_CZK_OK All is well * @retval C_CZK_BADARGS Invalid parameters were supplied */ -C_KZG_RET expand_root_of_unity(fr_t *out, const fr_t *root, uint64_t width) { +static C_KZG_RET expand_root_of_unity(fr_t *out, const fr_t *root, uint64_t width) { out[0] = fr_one; out[1] = *root; @@ -99,3 +99,81 @@ void free_fft_settings(FFTSettings *fs) { free(fs->reverse_roots_of_unity); fs->max_width = 0; } + +#ifdef KZGTEST + +#include "../inc/acutest.h" +#include "test_util.h" +#include "utility.h" + +#define NUM_ROOTS 32 + +void roots_of_unity_is_the_expected_size(void) { + TEST_CHECK(NUM_ROOTS == sizeof(scale2_root_of_unity) / sizeof(scale2_root_of_unity[0])); +} + +void roots_of_unity_out_of_bounds_fails(void) { + FFTSettings fs; + TEST_CHECK(C_KZG_BADARGS == new_fft_settings(&fs, NUM_ROOTS)); +} + +void roots_of_unity_are_plausible(void) { + fr_t r; + for (int i = 0; i < NUM_ROOTS; i++) { + fr_from_uint64s(&r, scale2_root_of_unity[i]); + for (int j = 0; j < i; j++) { + fr_sqr(&r, &r); + } + TEST_CHECK(true == fr_is_one(&r)); + TEST_MSG("Root %d failed", i); + } +} + +void expand_roots_is_plausible(void) { + // Just test one (largeish) value of scale + unsigned int scale = 15; + unsigned int width = 1 << scale; + fr_t root, expanded[width + 1], prod; + + // Initialise + fr_from_uint64s(&root, scale2_root_of_unity[scale]); + TEST_CHECK(expand_root_of_unity(expanded, &root, width) == C_KZG_OK); + + // Verify - each pair should multiply to one + TEST_CHECK(true == fr_is_one(expanded + 0)); + TEST_CHECK(true == fr_is_one(expanded + width)); + for (unsigned int i = 1; i <= width / 2; i++) { + fr_mul(&prod, expanded + i, expanded + width - i); + TEST_CHECK(true == fr_is_one(&prod)); + } +} + +void new_fft_settings_is_plausible(void) { + // Just test one (largeish) value of scale + int scale = 21; + unsigned int width = 1 << scale; + fr_t prod; + FFTSettings s; + + TEST_CHECK(new_fft_settings(&s, scale) == C_KZG_OK); + + // Verify - each pair should multiply to one + for (unsigned int i = 1; i <= width; i++) { + fr_mul(&prod, s.expanded_roots_of_unity + i, s.reverse_roots_of_unity + i); + TEST_CHECK(true == fr_is_one(&prod)); + } + + free_fft_settings(&s); +} + +TEST_LIST = { + {"FFT_COMMON_TEST", title}, + {"roots_of_unity_is_the_expected_size", roots_of_unity_is_the_expected_size}, + {"roots_of_unity_out_of_bounds_fails", roots_of_unity_out_of_bounds_fails}, + {"roots_of_unity_are_plausible", roots_of_unity_are_plausible}, + {"expand_roots_is_plausible", expand_roots_is_plausible}, + {"new_fft_settings_is_plausible", new_fft_settings_is_plausible}, + {NULL, NULL} /* zero record marks the end of the list */ +}; + +#endif // KZGTEST \ No newline at end of file diff --git a/src/fft_common.h b/src/fft_common.h index 90a1425..6e933fe 100644 --- a/src/fft_common.h +++ b/src/fft_common.h @@ -87,7 +87,6 @@ typedef struct { fr_t *reverse_roots_of_unity; /**< Descending powers of the root of unity, size `width + 1`. */ } FFTSettings; -C_KZG_RET expand_root_of_unity(fr_t *out, const fr_t *root, uint64_t width); C_KZG_RET new_fft_settings(FFTSettings *s, unsigned int max_scale); void free_fft_settings(FFTSettings *s); diff --git a/src/fft_common_test.c b/src/fft_common_test.c deleted file mode 100644 index d77effd..0000000 --- a/src/fft_common_test.c +++ /dev/null @@ -1,90 +0,0 @@ -/* - * Copyright 2021 Benjamin Edgington - * - * Licensed under the Apache License, Version 2.0 (the "License"); - * you may not use this file except in compliance with the License. - * You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -#include "../inc/acutest.h" -#include "test_util.h" -#include "fft_common.h" -#include "utility.h" - -#define NUM_ROOTS 32 - -void roots_of_unity_is_the_expected_size(void) { - TEST_CHECK(NUM_ROOTS == sizeof(scale2_root_of_unity) / sizeof(scale2_root_of_unity[0])); -} - -void roots_of_unity_out_of_bounds_fails(void) { - FFTSettings fs; - TEST_CHECK(C_KZG_BADARGS == new_fft_settings(&fs, NUM_ROOTS)); -} - -void roots_of_unity_are_plausible(void) { - fr_t r; - for (int i = 0; i < NUM_ROOTS; i++) { - fr_from_uint64s(&r, scale2_root_of_unity[i]); - for (int j = 0; j < i; j++) { - fr_sqr(&r, &r); - } - TEST_CHECK(true == fr_is_one(&r)); - TEST_MSG("Root %d failed", i); - } -} - -void expand_roots_is_plausible(void) { - // Just test one (largeish) value of scale - unsigned int scale = 15; - unsigned int width = 1 << scale; - fr_t root, expanded[width + 1], prod; - - // Initialise - fr_from_uint64s(&root, scale2_root_of_unity[scale]); - TEST_CHECK(expand_root_of_unity(expanded, &root, width) == C_KZG_OK); - - // Verify - each pair should multiply to one - TEST_CHECK(true == fr_is_one(expanded + 0)); - TEST_CHECK(true == fr_is_one(expanded + width)); - for (unsigned int i = 1; i <= width / 2; i++) { - fr_mul(&prod, expanded + i, expanded + width - i); - TEST_CHECK(true == fr_is_one(&prod)); - } -} - -void new_fft_settings_is_plausible(void) { - // Just test one (largeish) value of scale - int scale = 21; - unsigned int width = 1 << scale; - fr_t prod; - FFTSettings s; - - TEST_CHECK(new_fft_settings(&s, scale) == C_KZG_OK); - - // Verify - each pair should multiply to one - for (unsigned int i = 1; i <= width; i++) { - fr_mul(&prod, s.expanded_roots_of_unity + i, s.reverse_roots_of_unity + i); - TEST_CHECK(true == fr_is_one(&prod)); - } - - free_fft_settings(&s); -} - -TEST_LIST = { - {"FFT_COMMON_TEST", title}, - {"roots_of_unity_is_the_expected_size", roots_of_unity_is_the_expected_size}, - {"roots_of_unity_out_of_bounds_fails", roots_of_unity_out_of_bounds_fails}, - {"roots_of_unity_are_plausible", roots_of_unity_are_plausible}, - {"expand_roots_is_plausible", expand_roots_is_plausible}, - {"new_fft_settings_is_plausible", new_fft_settings_is_plausible}, - {NULL, NULL} /* zero record marks the end of the list */ -}; diff --git a/src/fft_fr.c b/src/fft_fr.c index 07f5f63..2db3d46 100644 --- a/src/fft_fr.c +++ b/src/fft_fr.c @@ -26,35 +26,8 @@ */ #include "fft_fr.h" -#include "c_kzg_util.h" #include "utility.h" -/** - * Slow Fourier Transform. - * - * This is simple, and ok for small sizes. It's mostly useful for testing. - * - * @param[out] out The results (array of length @p n) - * @param[in] in The input data (array of length @p n * @p stride) - * @param[in] stride The input data stride - * @param[in] roots Roots of unity (array of length @p n * @p roots_stride) - * @param[in] roots_stride The stride interval among the roots of unity - * @param[in] n Length of the FFT, must be a power of two - */ -void fft_fr_slow(fr_t *out, const fr_t *in, uint64_t stride, const fr_t *roots, uint64_t roots_stride, uint64_t n) { - fr_t v, last, jv, r; - for (uint64_t i = 0; i < n; i++) { - fr_mul(&last, &in[0], &roots[0]); - for (uint64_t j = 1; j < n; j++) { - jv = in[j * stride]; - r = roots[((i * j) % n) * roots_stride]; - fr_mul(&v, &jv, &r); - fr_add(&last, &last, &v); - } - out[i] = last; - } -} - /** * Fast Fourier Transform. * @@ -67,7 +40,8 @@ void fft_fr_slow(fr_t *out, const fr_t *in, uint64_t stride, const fr_t *roots, * @param[in] roots_stride The stride interval among the roots of unity * @param[in] n Length of the FFT, must be a power of two */ -void fft_fr_fast(fr_t *out, const fr_t *in, uint64_t stride, const fr_t *roots, uint64_t roots_stride, uint64_t n) { +static void fft_fr_fast(fr_t *out, const fr_t *in, uint64_t stride, const fr_t *roots, uint64_t roots_stride, + uint64_t n) { uint64_t half = n / 2; if (half > 0) { // Tunable parameter fft_fr_fast(out, in, stride * 2, roots, roots_stride * 2, half); @@ -111,3 +85,157 @@ C_KZG_RET fft_fr(fr_t *out, const fr_t *in, bool inverse, uint64_t n, const FFTS } return C_KZG_OK; } + +#ifdef KZGTEST + +#include "../inc/acutest.h" +#include "test_util.h" + +const uint64_t inv_fft_expected[][4] = { + {0x7fffffff80000008L, 0xa9ded2017fff2dffL, 0x199cec0404d0ec02L, 0x39f6d3a994cebea4L}, + {0xef296e7ffb8ca216L, 0xd5b902cbcef9c1b6L, 0xf06dfe5c7fca260dL, 0x13993b7d05187205L}, + {0xe930fdda2306c7d4L, 0x40e02aff48e2b16bL, 0x83a712d1dd818c8fL, 0x5dbc603bc53c7a3aL}, + {0xf9925986d0d25e90L, 0xcdf85d0a339d7782L, 0xee7a9a5f0410e423L, 0x2e0d216170831056L}, + {0x80007fff80000000L, 0x1fe05202bb00adffL, 0x6045d26b3fd26e6bL, 0x39f6d3a994cebea4L}, + {0x27325dd08ac4cee9L, 0xcbb94f168ddacca9L, 0x6843be68485784b1L, 0x5a6faf9039451673L}, + {0xe92ffdda2306c7d4L, 0x54dd2afcd2dfb16bL, 0xf6554603677e87beL, 0x5dbc603bc53c7a39L}, + {0x1cc772c9b57f126fL, 0xfb73f4d33d3116ddL, 0x4f9388c8d80abcf9L, 0x3ffbc9abcdda7821L}, + {0x7fffffff80000000L, 0xa9ded2017fff2dffL, 0x199cec0404d0ec02L, 0x39f6d3a994cebea4L}, + {0xe3388d354a80ed91L, 0x5849af2fc2cd4521L, 0xe3a64f3f31971b0bL, 0x33f1dda75bc30526L}, + {0x16d00224dcf9382cL, 0xfee079062d1eaa93L, 0x3ce49204a2235046L, 0x163147176461030eL}, + {0xd8cda22e753b3117L, 0x880454ec72238f55L, 0xcaf6199fc14a5353L, 0x197df7c2f05866d4L}, + {0x7fff7fff80000000L, 0x33dd520044fdadffL, 0xd2f4059cc9cf699aL, 0x39f6d3a994cebea3L}, + {0x066da6782f2da170L, 0x85c546f8cc60e47cL, 0x44bf3da90590f3e1L, 0x45e085f1b91a6cf1L}, + {0x16cf0224dcf9382cL, 0x12dd7903b71baa93L, 0xaf92c5362c204b76L, 0x163147176461030dL}, + {0x10d6917f04735deaL, 0x7e04a13731049a48L, 0x42cbd9ab89d7b1f7L, 0x60546bd624850b42L}}; + +/** + * Slow Fourier Transform. + * + * This is simple, and ok for small sizes. It's mostly useful for testing. + * + * @param[out] out The results (array of length @p n) + * @param[in] in The input data (array of length @p n * @p stride) + * @param[in] stride The input data stride + * @param[in] roots Roots of unity (array of length @p n * @p roots_stride) + * @param[in] roots_stride The stride interval among the roots of unity + * @param[in] n Length of the FFT, must be a power of two + */ +static void fft_fr_slow(fr_t *out, const fr_t *in, uint64_t stride, const fr_t *roots, uint64_t roots_stride, + uint64_t n) { + fr_t v, last, jv, r; + for (uint64_t i = 0; i < n; i++) { + fr_mul(&last, &in[0], &roots[0]); + for (uint64_t j = 1; j < n; j++) { + jv = in[j * stride]; + r = roots[((i * j) % n) * roots_stride]; + fr_mul(&v, &jv, &r); + fr_add(&last, &last, &v); + } + out[i] = last; + } +} + +void compare_sft_fft(void) { + // Initialise: ascending values of i (could be anything), and arbitrary size + unsigned int size = 12; + FFTSettings fs; + TEST_CHECK(new_fft_settings(&fs, size) == C_KZG_OK); + fr_t data[fs.max_width], out0[fs.max_width], out1[fs.max_width]; + for (int i = 0; i < fs.max_width; i++) { + fr_from_uint64(data + i, i); + } + + // Do both fast and slow transforms + fft_fr_slow(out0, data, 1, fs.expanded_roots_of_unity, 1, fs.max_width); + fft_fr_fast(out1, data, 1, fs.expanded_roots_of_unity, 1, fs.max_width); + + // Verify the results are identical + for (int i = 0; i < fs.max_width; i++) { + TEST_CHECK(fr_equal(out0 + i, out1 + i)); + } + + free_fft_settings(&fs); +} + +void roundtrip_fft(void) { + // Initialise: ascending values of i, and arbitrary size + unsigned int size = 12; + FFTSettings fs; + TEST_CHECK(new_fft_settings(&fs, size) == C_KZG_OK); + fr_t data[fs.max_width], coeffs[fs.max_width]; + for (int i = 0; i < fs.max_width; i++) { + fr_from_uint64(data + i, i); + } + + // Forward and reverse FFT + TEST_CHECK(fft_fr(coeffs, data, false, fs.max_width, &fs) == C_KZG_OK); + TEST_CHECK(fft_fr(data, coeffs, true, fs.max_width, &fs) == C_KZG_OK); + + // Verify that the result is still ascending values of i + for (int i = 0; i < fs.max_width; i++) { + fr_t tmp; + fr_from_uint64(&tmp, i); + TEST_CHECK(fr_equal(&tmp, data + i)); + } + + free_fft_settings(&fs); +} + +void inverse_fft(void) { + // Initialise: ascending values of i + FFTSettings fs; + TEST_CHECK(new_fft_settings(&fs, 4) == C_KZG_OK); + fr_t data[fs.max_width], out[fs.max_width]; + for (int i = 0; i < fs.max_width; i++) { + fr_from_uint64(&data[i], i); + } + + // Inverst FFT + TEST_CHECK(fft_fr(out, data, true, fs.max_width, &fs) == C_KZG_OK); + + // Verify against the known result, `inv_fft_expected` + int n = sizeof inv_fft_expected / sizeof inv_fft_expected[0]; + TEST_CHECK(n == fs.max_width); + for (int i = 0; i < n; i++) { + fr_t expected; + fr_from_uint64s(&expected, inv_fft_expected[i]); + TEST_CHECK(fr_equal(&expected, &out[i])); + } + + free_fft_settings(&fs); +} + +void stride_fft(void) { + unsigned int size1 = 9, size2 = 12; + uint64_t width = size1 < size2 ? (uint64_t)1 << size1 : (uint64_t)1 << size2; + FFTSettings fs1, fs2; + TEST_CHECK(new_fft_settings(&fs1, size1) == C_KZG_OK); + TEST_CHECK(new_fft_settings(&fs2, size2) == C_KZG_OK); + fr_t data[width], coeffs1[width], coeffs2[width]; + for (int i = 0; i < width; i++) { + fr_from_uint64(data + i, i); + } + + TEST_CHECK(fft_fr(coeffs1, data, false, width, &fs1) == C_KZG_OK); + TEST_CHECK(fft_fr(coeffs2, data, false, width, &fs2) == C_KZG_OK); + + for (int i = 0; i < width; i++) { + TEST_CHECK(fr_equal(coeffs1 + i, coeffs2 + i)); + } + + free_fft_settings(&fs1); + free_fft_settings(&fs2); +} + +TEST_LIST = { + {"FFT_FR_TEST", title}, + {"compare_sft_fft", compare_sft_fft}, + {"roundtrip_fft", roundtrip_fft}, + {"inverse_fft", inverse_fft}, + {"stride_fft", stride_fft}, + + {NULL, NULL} /* zero record marks the end of the list */ +}; + +#endif // KZGTEST \ No newline at end of file diff --git a/src/fft_fr.h b/src/fft_fr.h index 90648f2..3a4ab3a 100644 --- a/src/fft_fr.h +++ b/src/fft_fr.h @@ -18,6 +18,4 @@ #include "fft_common.h" -void fft_fr_slow(fr_t *out, const fr_t *in, uint64_t stride, const fr_t *roots, uint64_t roots_stride, uint64_t n); -void fft_fr_fast(fr_t *out, const fr_t *in, uint64_t stride, const fr_t *roots, uint64_t roots_stride, uint64_t n); C_KZG_RET fft_fr(fr_t *out, const fr_t *in, bool inverse, uint64_t n, const FFTSettings *fs); diff --git a/src/fft_fr_test.c b/src/fft_fr_test.c deleted file mode 100644 index 4d03e55..0000000 --- a/src/fft_fr_test.c +++ /dev/null @@ -1,139 +0,0 @@ -/* - * Copyright 2021 Benjamin Edgington - * - * Licensed under the Apache License, Version 2.0 (the "License"); - * you may not use this file except in compliance with the License. - * You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -#include "../inc/acutest.h" -#include "test_util.h" -#include "fft_fr.h" - -const uint64_t inv_fft_expected[][4] = { - {0x7fffffff80000008L, 0xa9ded2017fff2dffL, 0x199cec0404d0ec02L, 0x39f6d3a994cebea4L}, - {0xef296e7ffb8ca216L, 0xd5b902cbcef9c1b6L, 0xf06dfe5c7fca260dL, 0x13993b7d05187205L}, - {0xe930fdda2306c7d4L, 0x40e02aff48e2b16bL, 0x83a712d1dd818c8fL, 0x5dbc603bc53c7a3aL}, - {0xf9925986d0d25e90L, 0xcdf85d0a339d7782L, 0xee7a9a5f0410e423L, 0x2e0d216170831056L}, - {0x80007fff80000000L, 0x1fe05202bb00adffL, 0x6045d26b3fd26e6bL, 0x39f6d3a994cebea4L}, - {0x27325dd08ac4cee9L, 0xcbb94f168ddacca9L, 0x6843be68485784b1L, 0x5a6faf9039451673L}, - {0xe92ffdda2306c7d4L, 0x54dd2afcd2dfb16bL, 0xf6554603677e87beL, 0x5dbc603bc53c7a39L}, - {0x1cc772c9b57f126fL, 0xfb73f4d33d3116ddL, 0x4f9388c8d80abcf9L, 0x3ffbc9abcdda7821L}, - {0x7fffffff80000000L, 0xa9ded2017fff2dffL, 0x199cec0404d0ec02L, 0x39f6d3a994cebea4L}, - {0xe3388d354a80ed91L, 0x5849af2fc2cd4521L, 0xe3a64f3f31971b0bL, 0x33f1dda75bc30526L}, - {0x16d00224dcf9382cL, 0xfee079062d1eaa93L, 0x3ce49204a2235046L, 0x163147176461030eL}, - {0xd8cda22e753b3117L, 0x880454ec72238f55L, 0xcaf6199fc14a5353L, 0x197df7c2f05866d4L}, - {0x7fff7fff80000000L, 0x33dd520044fdadffL, 0xd2f4059cc9cf699aL, 0x39f6d3a994cebea3L}, - {0x066da6782f2da170L, 0x85c546f8cc60e47cL, 0x44bf3da90590f3e1L, 0x45e085f1b91a6cf1L}, - {0x16cf0224dcf9382cL, 0x12dd7903b71baa93L, 0xaf92c5362c204b76L, 0x163147176461030dL}, - {0x10d6917f04735deaL, 0x7e04a13731049a48L, 0x42cbd9ab89d7b1f7L, 0x60546bd624850b42L}}; - -void compare_sft_fft(void) { - // Initialise: ascending values of i (could be anything), and arbitrary size - unsigned int size = 12; - FFTSettings fs; - TEST_CHECK(new_fft_settings(&fs, size) == C_KZG_OK); - fr_t data[fs.max_width], out0[fs.max_width], out1[fs.max_width]; - for (int i = 0; i < fs.max_width; i++) { - fr_from_uint64(data + i, i); - } - - // Do both fast and slow transforms - fft_fr_slow(out0, data, 1, fs.expanded_roots_of_unity, 1, fs.max_width); - fft_fr_fast(out1, data, 1, fs.expanded_roots_of_unity, 1, fs.max_width); - - // Verify the results are identical - for (int i = 0; i < fs.max_width; i++) { - TEST_CHECK(fr_equal(out0 + i, out1 + i)); - } - - free_fft_settings(&fs); -} - -void roundtrip_fft(void) { - // Initialise: ascending values of i, and arbitrary size - unsigned int size = 12; - FFTSettings fs; - TEST_CHECK(new_fft_settings(&fs, size) == C_KZG_OK); - fr_t data[fs.max_width], coeffs[fs.max_width]; - for (int i = 0; i < fs.max_width; i++) { - fr_from_uint64(data + i, i); - } - - // Forward and reverse FFT - TEST_CHECK(fft_fr(coeffs, data, false, fs.max_width, &fs) == C_KZG_OK); - TEST_CHECK(fft_fr(data, coeffs, true, fs.max_width, &fs) == C_KZG_OK); - - // Verify that the result is still ascending values of i - for (int i = 0; i < fs.max_width; i++) { - fr_t tmp; - fr_from_uint64(&tmp, i); - TEST_CHECK(fr_equal(&tmp, data + i)); - } - - free_fft_settings(&fs); -} - -void inverse_fft(void) { - // Initialise: ascending values of i - FFTSettings fs; - TEST_CHECK(new_fft_settings(&fs, 4) == C_KZG_OK); - fr_t data[fs.max_width], out[fs.max_width]; - for (int i = 0; i < fs.max_width; i++) { - fr_from_uint64(&data[i], i); - } - - // Inverst FFT - TEST_CHECK(fft_fr(out, data, true, fs.max_width, &fs) == C_KZG_OK); - - // Verify against the known result, `inv_fft_expected` - int n = sizeof inv_fft_expected / sizeof inv_fft_expected[0]; - TEST_CHECK(n == fs.max_width); - for (int i = 0; i < n; i++) { - fr_t expected; - fr_from_uint64s(&expected, inv_fft_expected[i]); - TEST_CHECK(fr_equal(&expected, &out[i])); - } - - free_fft_settings(&fs); -} - -void stride_fft(void) { - unsigned int size1 = 9, size2 = 12; - uint64_t width = size1 < size2 ? (uint64_t)1 << size1 : (uint64_t)1 << size2; - FFTSettings fs1, fs2; - TEST_CHECK(new_fft_settings(&fs1, size1) == C_KZG_OK); - TEST_CHECK(new_fft_settings(&fs2, size2) == C_KZG_OK); - fr_t data[width], coeffs1[width], coeffs2[width]; - for (int i = 0; i < width; i++) { - fr_from_uint64(data + i, i); - } - - TEST_CHECK(fft_fr(coeffs1, data, false, width, &fs1) == C_KZG_OK); - TEST_CHECK(fft_fr(coeffs2, data, false, width, &fs2) == C_KZG_OK); - - for (int i = 0; i < width; i++) { - TEST_CHECK(fr_equal(coeffs1 + i, coeffs2 + i)); - } - - free_fft_settings(&fs1); - free_fft_settings(&fs2); -} - -TEST_LIST = { - {"FFT_FR_TEST", title}, - {"compare_sft_fft", compare_sft_fft}, - {"roundtrip_fft", roundtrip_fft}, - {"inverse_fft", inverse_fft}, - {"stride_fft", stride_fft}, - - {NULL, NULL} /* zero record marks the end of the list */ -}; diff --git a/src/fft_g1.c b/src/fft_g1.c index 2a38e94..99c87c6 100644 --- a/src/fft_g1.c +++ b/src/fft_g1.c @@ -26,36 +26,8 @@ */ #include "fft_g1.h" -#include "c_kzg_util.h" #include "utility.h" -/** - * Slow Fourier Transform. - * - * This is simple, and ok for small sizes. It's mostly useful for testing. - * - * @param[out] out The results (array of length @p n) - * @param[in] in The input data (array of length @p n * @p stride) - * @param[in] stride The input data stride - * @param[in] roots Roots of unity (array of length @p n * @p roots_stride) - * @param[in] roots_stride The stride interval among the roots of unity - * @param[in] n Length of the FFT, must be a power of two - */ -void fft_g1_slow(g1_t *out, const g1_t *in, uint64_t stride, const fr_t *roots, uint64_t roots_stride, uint64_t n) { - g1_t v, last, jv; - fr_t r; - for (uint64_t i = 0; i < n; i++) { - g1_mul(&last, &in[0], &roots[0]); - for (uint64_t j = 1; j < n; j++) { - jv = in[j * stride]; - r = roots[((i * j) % n) * roots_stride]; - g1_mul(&v, &jv, &r); - g1_add_or_dbl(&last, &last, &v); - } - out[i] = last; - } -} - /** * Fast Fourier Transform. * @@ -68,7 +40,8 @@ void fft_g1_slow(g1_t *out, const g1_t *in, uint64_t stride, const fr_t *roots, * @param[in] roots_stride The stride interval among the roots of unity * @param[in] n Length of the FFT, must be a power of two */ -void fft_g1_fast(g1_t *out, const g1_t *in, uint64_t stride, const fr_t *roots, uint64_t roots_stride, uint64_t n) { +static void fft_g1_fast(g1_t *out, const g1_t *in, uint64_t stride, const fr_t *roots, uint64_t roots_stride, + uint64_t n) { uint64_t half = n / 2; if (half > 0) { // Tunable parameter fft_g1_fast(out, in, stride * 2, roots, roots_stride * 2, half); @@ -112,3 +85,117 @@ C_KZG_RET fft_g1(g1_t *out, const g1_t *in, bool inverse, uint64_t n, const FFTS } return C_KZG_OK; } + +#ifdef KZGTEST + +#include "../inc/acutest.h" +#include "test_util.h" +#include "fft_g1.h" + +static void make_data(g1_t *out, uint64_t n) { + // Multiples of g1_gen + if (n == 0) return; + *out = g1_generator; + for (int i = 1; i < n; i++) { + g1_add_or_dbl(out + i, out + i - 1, &g1_generator); + } +} + +/** + * Slow Fourier Transform. + * + * This is simple, and ok for small sizes. It's mostly useful for testing. + * + * @param[out] out The results (array of length @p n) + * @param[in] in The input data (array of length @p n * @p stride) + * @param[in] stride The input data stride + * @param[in] roots Roots of unity (array of length @p n * @p roots_stride) + * @param[in] roots_stride The stride interval among the roots of unity + * @param[in] n Length of the FFT, must be a power of two + */ +static void fft_g1_slow(g1_t *out, const g1_t *in, uint64_t stride, const fr_t *roots, uint64_t roots_stride, + uint64_t n) { + g1_t v, last, jv; + fr_t r; + for (uint64_t i = 0; i < n; i++) { + g1_mul(&last, &in[0], &roots[0]); + for (uint64_t j = 1; j < n; j++) { + jv = in[j * stride]; + r = roots[((i * j) % n) * roots_stride]; + g1_mul(&v, &jv, &r); + g1_add_or_dbl(&last, &last, &v); + } + out[i] = last; + } +} + +void compare_sft_fft(void) { + // Initialise: arbitrary size + unsigned int size = 6; + FFTSettings fs; + TEST_CHECK(new_fft_settings(&fs, size) == C_KZG_OK); + g1_t data[fs.max_width], slow[fs.max_width], fast[fs.max_width]; + make_data(data, fs.max_width); + + // Do both fast and slow transforms + fft_g1_slow(slow, data, 1, fs.expanded_roots_of_unity, 1, fs.max_width); + fft_g1_fast(fast, data, 1, fs.expanded_roots_of_unity, 1, fs.max_width); + + // Verify the results are identical + for (int i = 0; i < fs.max_width; i++) { + TEST_CHECK(g1_equal(slow + i, fast + i)); + } + + free_fft_settings(&fs); +} + +void roundtrip_fft(void) { + // Initialise: arbitrary size + unsigned int size = 10; + FFTSettings fs; + TEST_CHECK(new_fft_settings(&fs, size) == C_KZG_OK); + g1_t expected[fs.max_width], data[fs.max_width], coeffs[fs.max_width]; + make_data(expected, fs.max_width); + make_data(data, fs.max_width); + + // Forward and reverse FFT + TEST_CHECK(fft_g1(coeffs, data, false, fs.max_width, &fs) == C_KZG_OK); + TEST_CHECK(fft_g1(data, coeffs, true, fs.max_width, &fs) == C_KZG_OK); + + // Verify that the result is still ascending values of i + for (int i = 0; i < fs.max_width; i++) { + TEST_CHECK(g1_equal(expected + i, data + i)); + } + + free_fft_settings(&fs); +} + +void stride_fft(void) { + unsigned int size1 = 9, size2 = 12; + uint64_t width = size1 < size2 ? (uint64_t)1 << size1 : (uint64_t)1 << size2; + FFTSettings fs1, fs2; + TEST_CHECK(new_fft_settings(&fs1, size1) == C_KZG_OK); + TEST_CHECK(new_fft_settings(&fs2, size2) == C_KZG_OK); + g1_t data[width], coeffs1[width], coeffs2[width]; + make_data(data, width); + + TEST_CHECK(fft_g1(coeffs1, data, false, width, &fs1) == C_KZG_OK); + TEST_CHECK(fft_g1(coeffs2, data, false, width, &fs2) == C_KZG_OK); + + for (int i = 0; i < width; i++) { + TEST_CHECK(g1_equal(coeffs1 + i, coeffs2 + i)); + } + + free_fft_settings(&fs1); + free_fft_settings(&fs2); +} + +TEST_LIST = { + {"FFT_G1_TEST", title}, + {"compare_sft_fft", compare_sft_fft}, + {"roundtrip_fft", roundtrip_fft}, + {"stride_fft", stride_fft}, + {NULL, NULL} /* zero record marks the end of the list */ +}; + +#endif // KZGTEST \ No newline at end of file diff --git a/src/fft_g1.h b/src/fft_g1.h index 1e34654..e0c7bd3 100644 --- a/src/fft_g1.h +++ b/src/fft_g1.h @@ -18,6 +18,4 @@ #include "fft_common.h" -void fft_g1_slow(g1_t *out, const g1_t *in, uint64_t stride, const fr_t *roots, uint64_t roots_stride, uint64_t n); -void fft_g1_fast(g1_t *out, const g1_t *in, uint64_t stride, const fr_t *roots, uint64_t roots_stride, uint64_t n); C_KZG_RET fft_g1(g1_t *out, const g1_t *in, bool inverse, uint64_t n, const FFTSettings *fs); diff --git a/src/fft_g1_test.c b/src/fft_g1_test.c deleted file mode 100644 index 7c812e1..0000000 --- a/src/fft_g1_test.c +++ /dev/null @@ -1,97 +0,0 @@ -/* - * Copyright 2021 Benjamin Edgington - * - * Licensed under the Apache License, Version 2.0 (the "License"); - * you may not use this file except in compliance with the License. - * You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -#include "../inc/acutest.h" -#include "test_util.h" -#include "fft_g1.h" - -void make_data(g1_t *out, uint64_t n) { - // Multiples of g1_gen - if (n == 0) return; - *out = g1_generator; - for (int i = 1; i < n; i++) { - g1_add_or_dbl(out + i, out + i - 1, &g1_generator); - } -} - -void compare_sft_fft(void) { - // Initialise: arbitrary size - unsigned int size = 6; - FFTSettings fs; - TEST_CHECK(new_fft_settings(&fs, size) == C_KZG_OK); - g1_t data[fs.max_width], slow[fs.max_width], fast[fs.max_width]; - make_data(data, fs.max_width); - - // Do both fast and slow transforms - fft_g1_slow(slow, data, 1, fs.expanded_roots_of_unity, 1, fs.max_width); - fft_g1_fast(fast, data, 1, fs.expanded_roots_of_unity, 1, fs.max_width); - - // Verify the results are identical - for (int i = 0; i < fs.max_width; i++) { - TEST_CHECK(g1_equal(slow + i, fast + i)); - } - - free_fft_settings(&fs); -} - -void roundtrip_fft(void) { - // Initialise: arbitrary size - unsigned int size = 10; - FFTSettings fs; - TEST_CHECK(new_fft_settings(&fs, size) == C_KZG_OK); - g1_t expected[fs.max_width], data[fs.max_width], coeffs[fs.max_width]; - make_data(expected, fs.max_width); - make_data(data, fs.max_width); - - // Forward and reverse FFT - TEST_CHECK(fft_g1(coeffs, data, false, fs.max_width, &fs) == C_KZG_OK); - TEST_CHECK(fft_g1(data, coeffs, true, fs.max_width, &fs) == C_KZG_OK); - - // Verify that the result is still ascending values of i - for (int i = 0; i < fs.max_width; i++) { - TEST_CHECK(g1_equal(expected + i, data + i)); - } - - free_fft_settings(&fs); -} - -void stride_fft(void) { - unsigned int size1 = 9, size2 = 12; - uint64_t width = size1 < size2 ? (uint64_t)1 << size1 : (uint64_t)1 << size2; - FFTSettings fs1, fs2; - TEST_CHECK(new_fft_settings(&fs1, size1) == C_KZG_OK); - TEST_CHECK(new_fft_settings(&fs2, size2) == C_KZG_OK); - g1_t data[width], coeffs1[width], coeffs2[width]; - make_data(data, width); - - TEST_CHECK(fft_g1(coeffs1, data, false, width, &fs1) == C_KZG_OK); - TEST_CHECK(fft_g1(coeffs2, data, false, width, &fs2) == C_KZG_OK); - - for (int i = 0; i < width; i++) { - TEST_CHECK(g1_equal(coeffs1 + i, coeffs2 + i)); - } - - free_fft_settings(&fs1); - free_fft_settings(&fs2); -} - -TEST_LIST = { - {"FFT_G1_TEST", title}, - {"compare_sft_fft", compare_sft_fft}, - {"roundtrip_fft", roundtrip_fft}, - {"stride_fft", stride_fft}, - {NULL, NULL} /* zero record marks the end of the list */ -}; diff --git a/src/fk20_proofs.c b/src/fk20_proofs.c index e2478fe..b56c0d6 100644 --- a/src/fk20_proofs.c +++ b/src/fk20_proofs.c @@ -40,7 +40,7 @@ * @retval C_CZK_ERROR An internal error occurred * @retval C_CZK_MALLOC Memory allocation failed */ -C_KZG_RET toeplitz_part_1(g1_t *out, const g1_t *x, uint64_t n, const FFTSettings *fs) { +static C_KZG_RET toeplitz_part_1(g1_t *out, const g1_t *x, uint64_t n, const FFTSettings *fs) { uint64_t n2 = n * 2; g1_t *x_ext; @@ -70,7 +70,7 @@ C_KZG_RET toeplitz_part_1(g1_t *out, const g1_t *x, uint64_t n, const FFTSetting * @retval C_CZK_ERROR An internal error occurred * @retval C_CZK_MALLOC Memory allocation failed */ -C_KZG_RET toeplitz_part_2(g1_t *out, const poly *toeplitz_coeffs, const g1_t *x_ext_fft, const FFTSettings *fs) { +static C_KZG_RET toeplitz_part_2(g1_t *out, const poly *toeplitz_coeffs, const g1_t *x_ext_fft, const FFTSettings *fs) { fr_t *toeplitz_coeffs_fft; // CHECK(toeplitz_coeffs->length == fk->x_ext_fft_len); // TODO: how to implement? @@ -96,7 +96,7 @@ C_KZG_RET toeplitz_part_2(g1_t *out, const poly *toeplitz_coeffs, const g1_t *x_ * @retval C_CZK_OK All is well * @retval C_CZK_ERROR An internal error occurred */ -C_KZG_RET toeplitz_part_3(g1_t *out, const g1_t *h_ext_fft, uint64_t n2, const FFTSettings *fs) { +static C_KZG_RET toeplitz_part_3(g1_t *out, const g1_t *h_ext_fft, uint64_t n2, const FFTSettings *fs) { uint64_t n = n2 / 2; TRY(fft_g1(out, h_ext_fft, true, n2, fs)); @@ -122,7 +122,7 @@ C_KZG_RET toeplitz_part_3(g1_t *out, const g1_t *h_ext_fft, uint64_t n2, const F * @retval C_CZK_BADARGS Invalid parameters were supplied * @retval C_CZK_MALLOC Memory allocation failed */ -C_KZG_RET toeplitz_coeffs_stride(poly *out, const poly *in, uint64_t offset, uint64_t stride) { +static C_KZG_RET toeplitz_coeffs_stride(poly *out, const poly *in, uint64_t offset, uint64_t stride) { uint64_t n = in->length, k, k2; CHECK(stride > 0); @@ -151,7 +151,7 @@ C_KZG_RET toeplitz_coeffs_stride(poly *out, const poly *in, uint64_t offset, uin * @retval C_CZK_OK All is well * @retval C_CZK_MALLOC Memory allocation failed */ -C_KZG_RET toeplitz_coeffs_step(poly *out, const poly *in) { +static C_KZG_RET toeplitz_coeffs_step(poly *out, const poly *in) { return toeplitz_coeffs_stride(out, in, 0, 1); } @@ -174,7 +174,7 @@ C_KZG_RET toeplitz_coeffs_step(poly *out, const poly *in) { * @retval C_CZK_ERROR An internal error occurred * @retval C_CZK_MALLOC Memory allocation failed */ -C_KZG_RET fk20_single_da_opt(g1_t *out, const poly *p, const FK20SingleSettings *fk) { +static C_KZG_RET fk20_single_da_opt(g1_t *out, const poly *p, const FK20SingleSettings *fk) { uint64_t n = p->length, n2 = n * 2; g1_t *h, *h_ext_fft; poly toeplitz_coeffs; @@ -291,7 +291,7 @@ C_KZG_RET fk20_compute_proof_multi(g1_t *out, const poly *p, const FK20MultiSett * @param[in] p The polynomial, length `n` * @param[in] fk FK20 multi settings previously initialised by #new_fk20_multi_settings */ -C_KZG_RET fk20_multi_da_opt(g1_t *out, const poly *p, const FK20MultiSettings *fk) { +static C_KZG_RET fk20_multi_da_opt(g1_t *out, const poly *p, const FK20MultiSettings *fk) { uint64_t n = p->length, n2 = n * 2, k, k2; g1_t *h_ext_fft, *h_ext_fft_file, *h; poly toeplitz_coeffs; @@ -467,3 +467,279 @@ void free_fk20_multi_settings(FK20MultiSettings *fk) { fk->chunk_len = 0; fk->length = 0; } + +#ifdef KZGTEST + +#include "../inc/acutest.h" +#include "test_util.h" + +void fk_single(void) { + // Our polynomial: degree 15, 16 coefficients + uint64_t coeffs[] = {1, 2, 3, 4, 7, 7, 7, 7, 13, 13, 13, 13, 13, 13, 13, 13}; + int poly_len = sizeof coeffs / sizeof coeffs[0]; + + // The FFT settings size + uint64_t n = 5, n_len = (uint64_t)1 << n; + + FFTSettings fs; + KZGSettings ks; + FK20SingleSettings fk; + uint64_t secrets_len = n_len + 1; + g1_t s1[secrets_len]; + g2_t s2[secrets_len]; + poly p; + g1_t commitment, all_proofs[2 * poly_len], proof; + fr_t x, y; + bool result; + + TEST_CHECK(n_len >= 2 * poly_len); + TEST_CHECK(new_poly(&p, poly_len) == C_KZG_OK); + for (uint64_t i = 0; i < poly_len; i++) { + fr_from_uint64(&p.coeffs[i], coeffs[i]); + } + + // Initialise the secrets and data structures + generate_trusted_setup(s1, s2, &secret, secrets_len); + TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, n)); + TEST_CHECK(C_KZG_OK == new_kzg_settings(&ks, s1, s2, secrets_len, &fs)); + TEST_CHECK(C_KZG_OK == new_fk20_single_settings(&fk, 2 * poly_len, &ks)); + + // Commit to the polynomial + TEST_CHECK(C_KZG_OK == commit_to_poly(&commitment, &p, &ks)); + + // 1. First with `da_using_fk20_single` + + // Generate the proofs + TEST_CHECK(da_using_fk20_single(all_proofs, &p, &fk) == C_KZG_OK); + + // Verify the proof at each root of unity + for (uint64_t i = 0; i < 2 * poly_len; i++) { + x = fs.expanded_roots_of_unity[i]; + eval_poly(&y, &p, &x); + proof = all_proofs[reverse_bits_limited(2 * poly_len, i)]; + + TEST_CHECK(C_KZG_OK == check_proof_single(&result, &commitment, &proof, &x, &y, &ks)); + TEST_CHECK(true == result); + } + + // 2. Exactly the same thing again with `fk20_single_da_opt` + + // Generate the proofs + TEST_CHECK(fk20_single_da_opt(all_proofs, &p, &fk) == C_KZG_OK); + + // Verify the proof at each root of unity + for (uint64_t i = 0; i < 2 * poly_len; i++) { + x = fs.expanded_roots_of_unity[i]; + eval_poly(&y, &p, &x); + proof = all_proofs[i]; + + TEST_CHECK(C_KZG_OK == check_proof_single(&result, &commitment, &proof, &x, &y, &ks)); + TEST_CHECK(true == result); + } + + // Clean up + free_poly(&p); + free_fft_settings(&fs); + free_kzg_settings(&ks); + free_fk20_single_settings(&fk); +} + +void fk_single_strided(void) { + // Our polynomial: degree 15, 16 coefficients + uint64_t coeffs[] = {1, 2, 3, 4, 7, 7, 7, 7, 13, 13, 13, 13, 13, 13, 13, 13}; + int poly_len = sizeof coeffs / sizeof coeffs[0]; + + // We can set up the FFTs for bigger widths if we wish. + // This is a useful canary for issues elsewhere in the code. + uint64_t n = 8, n_len = (uint64_t)1 << n; + uint64_t stride = n_len / (2 * poly_len); + + FFTSettings fs; + KZGSettings ks; + FK20SingleSettings fk; + uint64_t secrets_len = n_len + 1; + g1_t s1[secrets_len]; + g2_t s2[secrets_len]; + poly p; + g1_t commitment, all_proofs[2 * poly_len], proof; + fr_t x, y; + bool result; + + TEST_CHECK(n_len >= 2 * poly_len); + TEST_CHECK(new_poly(&p, poly_len) == C_KZG_OK); + for (uint64_t i = 0; i < poly_len; i++) { + fr_from_uint64(&p.coeffs[i], coeffs[i]); + } + + // Initialise the secrets and data structures + generate_trusted_setup(s1, s2, &secret, secrets_len); + TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, n)); + TEST_CHECK(C_KZG_OK == new_kzg_settings(&ks, s1, s2, secrets_len, &fs)); + TEST_CHECK(C_KZG_OK == new_fk20_single_settings(&fk, 2 * poly_len, &ks)); + + // Commit to the polynomial + TEST_CHECK(C_KZG_OK == commit_to_poly(&commitment, &p, &ks)); + + // Generate the proofs + TEST_CHECK(da_using_fk20_single(all_proofs, &p, &fk) == C_KZG_OK); + + // Verify the proof at each root of unity + for (uint64_t i = 0; i < 2 * poly_len; i++) { + x = fs.expanded_roots_of_unity[i * stride]; + eval_poly(&y, &p, &x); + proof = all_proofs[reverse_bits_limited(2 * poly_len, i)]; + + TEST_CHECK(C_KZG_OK == check_proof_single(&result, &commitment, &proof, &x, &y, &ks)); + TEST_CHECK(true == result); + } + + // Clean up + free_poly(&p); + free_fft_settings(&fs); + free_kzg_settings(&ks); + free_fk20_single_settings(&fk); +} + +void fk_multi_settings(void) { + FFTSettings fs; + KZGSettings ks; + FK20MultiSettings fk; + uint64_t n = 5; + uint64_t secrets_len = 33; + g1_t s1[secrets_len]; + g2_t s2[secrets_len]; + + // Initialise the secrets and data structures + generate_trusted_setup(s1, s2, &secret, secrets_len); + TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, n)); + TEST_CHECK(C_KZG_OK == new_kzg_settings(&ks, s1, s2, secrets_len, &fs)); + TEST_CHECK(C_KZG_OK == new_fk20_multi_settings(&fk, 32, 4, &ks)); + + // Don't do anything. Run this with `valgrind` to check that memory is correctly allocated and freed. + + free_fft_settings(&fs); + free_kzg_settings(&ks); + free_fk20_multi_settings(&fk); +} + +void fk_multi_0(void) { + FFTSettings fs; + KZGSettings ks; + FK20MultiSettings fk; + uint64_t n, chunk_len, chunk_count; + uint64_t secrets_len; + g1_t *s1; + g2_t *s2; + poly p; + uint64_t vv[] = {1, 2, 3, 4, 7, 8, 9, 10, 13, 14, 1, 15, 1, 1000, 134, 33}; + g1_t commitment; + g1_t *all_proofs; + fr_t *extended_coeffs, *extended_coeffs_fft; + fr_t *ys, *ys2; + uint64_t domain_stride; + + chunk_len = 16; + chunk_count = 32; + n = chunk_len * chunk_count; + secrets_len = 2 * n; + + TEST_CHECK(C_KZG_OK == new_g1_array(&s1, secrets_len)); + TEST_CHECK(C_KZG_OK == new_g2_array(&s2, secrets_len)); + + generate_trusted_setup(s1, s2, &secret, secrets_len); + TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 4 + 5 + 1)); + TEST_CHECK(C_KZG_OK == new_kzg_settings(&ks, s1, s2, secrets_len, &fs)); + TEST_CHECK(C_KZG_OK == new_fk20_multi_settings(&fk, n * 2, chunk_len, &ks)); + + // Create a test polynomial: 512 coefficients + TEST_CHECK(C_KZG_OK == new_poly(&p, n)); + for (int i = 0; i < chunk_count; i++) { + for (int j = 0; j < chunk_len; j++) { + uint64_t v = vv[j]; + if (j == 3) v += i; + if (j == 5) v += i * i; + fr_from_uint64(&p.coeffs[i * chunk_len + j], v); + } + fr_negate(&p.coeffs[i * chunk_len + 12], &p.coeffs[i * chunk_len + 12]); + fr_negate(&p.coeffs[i * chunk_len + 14], &p.coeffs[i * chunk_len + 14]); + } + + TEST_CHECK(C_KZG_OK == commit_to_poly(&commitment, &p, &ks)); + + // Compute the multi proofs, assuming that the polynomial will be extended with zeros + TEST_CHECK(C_KZG_OK == new_g1_array(&all_proofs, 2 * chunk_count)); + TEST_CHECK(C_KZG_OK == da_using_fk20_multi(all_proofs, &p, &fk)); + + // Now actually extend the polynomial with zeros + TEST_CHECK(C_KZG_OK == new_fr_array(&extended_coeffs, 2 * n)); + for (uint64_t i = 0; i < n; i++) { + extended_coeffs[i] = p.coeffs[i]; + } + for (uint64_t i = n; i < 2 * n; i++) { + extended_coeffs[i] = fr_zero; + } + TEST_CHECK(C_KZG_OK == new_fr_array(&extended_coeffs_fft, 2 * n)); + TEST_CHECK(C_KZG_OK == fft_fr(extended_coeffs_fft, extended_coeffs, false, 2 * n, &fs)); + TEST_CHECK(C_KZG_OK == reverse_bit_order(extended_coeffs_fft, sizeof extended_coeffs_fft[0], 2 * n)); + + // Verify the proofs + TEST_CHECK(C_KZG_OK == new_fr_array(&ys, chunk_len)); + TEST_CHECK(C_KZG_OK == new_fr_array(&ys2, chunk_len)); + domain_stride = fs.max_width / (2 * n); + for (uint64_t pos = 0; pos < 2 * chunk_count; pos++) { + uint64_t domain_pos, stride; + fr_t x; + bool result; + + domain_pos = reverse_bits_limited(2 * chunk_count, pos); + x = fs.expanded_roots_of_unity[domain_pos * domain_stride]; + + // The ys from the extended coeffients + for (uint64_t i = 0; i < chunk_len; i++) { + ys[i] = extended_coeffs_fft[chunk_len * pos + i]; + } + TEST_CHECK(C_KZG_OK == reverse_bit_order(ys, sizeof ys[0], chunk_len)); + + // Now recreate the ys by evaluating the polynomial in the sub-domain range + stride = fs.max_width / chunk_len; + for (uint64_t i = 0; i < chunk_len; i++) { + fr_t z; + fr_mul(&z, &x, &fs.expanded_roots_of_unity[i * stride]); + eval_poly(&ys2[i], &p, &z); + } + + // ys and ys2 should be equal + for (uint64_t i = 0; i < chunk_len; i++) { + TEST_CHECK(fr_equal(&ys[i], &ys2[i])); + } + + // Verify this proof + TEST_CHECK(C_KZG_OK == check_proof_multi(&result, &commitment, &all_proofs[pos], &x, ys, chunk_len, &ks)); + TEST_CHECK(true == result); + } + + free_poly(&p); + free(all_proofs); + free(extended_coeffs); + free(extended_coeffs_fft); + free(ys); + free(ys2); + free(s1); + free(s2); + free_fft_settings(&fs); + free_kzg_settings(&ks); + free_fk20_multi_settings(&fk); +} + +// TODO: compare results of fk20_multi_da_opt() and fk20_compute_proof_multi() + +TEST_LIST = { + {"FK20_PROOFS_TEST", title}, + {"fk_single", fk_single}, + {"fk_single_strided", fk_single_strided}, + {"fk_multi_settings", fk_multi_settings}, + {"fk_multi_0", fk_multi_0}, + {NULL, NULL} /* zero record marks the end of the list */ +}; + +#endif \ No newline at end of file diff --git a/src/fk20_proofs.h b/src/fk20_proofs.h index 8da72ad..540f7f0 100644 --- a/src/fk20_proofs.h +++ b/src/fk20_proofs.h @@ -39,14 +39,7 @@ typedef struct { uint64_t length; /**< TODO */ } FK20MultiSettings; -C_KZG_RET toeplitz_part_1(g1_t *out, const g1_t *x, uint64_t n, const FFTSettings *fs); -C_KZG_RET toeplitz_part_2(g1_t *out, const poly *toeplitz_coeffs, const g1_t *x_ext_fft, const FFTSettings *fs); -C_KZG_RET toeplitz_part_3(g1_t *out, const g1_t *h_ext_fft, uint64_t n2, const FFTSettings *fs); -C_KZG_RET toeplitz_coeffs_stride(poly *out, const poly *in, uint64_t offset, uint64_t stride); -C_KZG_RET toeplitz_coeffs_step(poly *out, const poly *in); -C_KZG_RET fk20_single_da_opt(g1_t *out, const poly *p, const FK20SingleSettings *fk); C_KZG_RET da_using_fk20_single(g1_t *out, const poly *p, const FK20SingleSettings *fk); -C_KZG_RET fk20_multi_da_opt(g1_t *out, const poly *p, const FK20MultiSettings *fk); C_KZG_RET da_using_fk20_multi(g1_t *out, const poly *p, const FK20MultiSettings *fk); C_KZG_RET new_fk20_single_settings(FK20SingleSettings *fk, uint64_t n2, const KZGSettings *ks); C_KZG_RET new_fk20_multi_settings(FK20MultiSettings *fk, uint64_t n2, uint64_t chunk_len, const KZGSettings *ks); diff --git a/src/fk20_proofs_test.c b/src/fk20_proofs_test.c deleted file mode 100644 index b98c788..0000000 --- a/src/fk20_proofs_test.c +++ /dev/null @@ -1,290 +0,0 @@ -/* - * Copyright 2021 Benjamin Edgington - * - * Licensed under the Apache License, Version 2.0 (the "License"); - * you may not use this file except in compliance with the License. - * You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -#include "../inc/acutest.h" -#include "test_util.h" -#include "fk20_proofs.h" -#include "c_kzg_util.h" -#include "utility.h" - -void fk_single(void) { - // Our polynomial: degree 15, 16 coefficients - uint64_t coeffs[] = {1, 2, 3, 4, 7, 7, 7, 7, 13, 13, 13, 13, 13, 13, 13, 13}; - int poly_len = sizeof coeffs / sizeof coeffs[0]; - - // The FFT settings size - uint64_t n = 5, n_len = (uint64_t)1 << n; - - FFTSettings fs; - KZGSettings ks; - FK20SingleSettings fk; - uint64_t secrets_len = n_len + 1; - g1_t s1[secrets_len]; - g2_t s2[secrets_len]; - poly p; - g1_t commitment, all_proofs[2 * poly_len], proof; - fr_t x, y; - bool result; - - TEST_CHECK(n_len >= 2 * poly_len); - TEST_CHECK(new_poly(&p, poly_len) == C_KZG_OK); - for (uint64_t i = 0; i < poly_len; i++) { - fr_from_uint64(&p.coeffs[i], coeffs[i]); - } - - // Initialise the secrets and data structures - generate_trusted_setup(s1, s2, &secret, secrets_len); - TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, n)); - TEST_CHECK(C_KZG_OK == new_kzg_settings(&ks, s1, s2, secrets_len, &fs)); - TEST_CHECK(C_KZG_OK == new_fk20_single_settings(&fk, 2 * poly_len, &ks)); - - // Commit to the polynomial - TEST_CHECK(C_KZG_OK == commit_to_poly(&commitment, &p, &ks)); - - // 1. First with `da_using_fk20_single` - - // Generate the proofs - TEST_CHECK(da_using_fk20_single(all_proofs, &p, &fk) == C_KZG_OK); - - // Verify the proof at each root of unity - for (uint64_t i = 0; i < 2 * poly_len; i++) { - x = fs.expanded_roots_of_unity[i]; - eval_poly(&y, &p, &x); - proof = all_proofs[reverse_bits_limited(2 * poly_len, i)]; - - TEST_CHECK(C_KZG_OK == check_proof_single(&result, &commitment, &proof, &x, &y, &ks)); - TEST_CHECK(true == result); - } - - // 2. Exactly the same thing again with `fk20_single_da_opt` - - // Generate the proofs - TEST_CHECK(fk20_single_da_opt(all_proofs, &p, &fk) == C_KZG_OK); - - // Verify the proof at each root of unity - for (uint64_t i = 0; i < 2 * poly_len; i++) { - x = fs.expanded_roots_of_unity[i]; - eval_poly(&y, &p, &x); - proof = all_proofs[i]; - - TEST_CHECK(C_KZG_OK == check_proof_single(&result, &commitment, &proof, &x, &y, &ks)); - TEST_CHECK(true == result); - } - - // Clean up - free_poly(&p); - free_fft_settings(&fs); - free_kzg_settings(&ks); - free_fk20_single_settings(&fk); -} - -void fk_single_strided(void) { - // Our polynomial: degree 15, 16 coefficients - uint64_t coeffs[] = {1, 2, 3, 4, 7, 7, 7, 7, 13, 13, 13, 13, 13, 13, 13, 13}; - int poly_len = sizeof coeffs / sizeof coeffs[0]; - - // We can set up the FFTs for bigger widths if we wish. - // This is a useful canary for issues elsewhere in the code. - uint64_t n = 8, n_len = (uint64_t)1 << n; - uint64_t stride = n_len / (2 * poly_len); - - FFTSettings fs; - KZGSettings ks; - FK20SingleSettings fk; - uint64_t secrets_len = n_len + 1; - g1_t s1[secrets_len]; - g2_t s2[secrets_len]; - poly p; - g1_t commitment, all_proofs[2 * poly_len], proof; - fr_t x, y; - bool result; - - TEST_CHECK(n_len >= 2 * poly_len); - TEST_CHECK(new_poly(&p, poly_len) == C_KZG_OK); - for (uint64_t i = 0; i < poly_len; i++) { - fr_from_uint64(&p.coeffs[i], coeffs[i]); - } - - // Initialise the secrets and data structures - generate_trusted_setup(s1, s2, &secret, secrets_len); - TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, n)); - TEST_CHECK(C_KZG_OK == new_kzg_settings(&ks, s1, s2, secrets_len, &fs)); - TEST_CHECK(C_KZG_OK == new_fk20_single_settings(&fk, 2 * poly_len, &ks)); - - // Commit to the polynomial - TEST_CHECK(C_KZG_OK == commit_to_poly(&commitment, &p, &ks)); - - // Generate the proofs - TEST_CHECK(da_using_fk20_single(all_proofs, &p, &fk) == C_KZG_OK); - - // Verify the proof at each root of unity - for (uint64_t i = 0; i < 2 * poly_len; i++) { - x = fs.expanded_roots_of_unity[i * stride]; - eval_poly(&y, &p, &x); - proof = all_proofs[reverse_bits_limited(2 * poly_len, i)]; - - TEST_CHECK(C_KZG_OK == check_proof_single(&result, &commitment, &proof, &x, &y, &ks)); - TEST_CHECK(true == result); - } - - // Clean up - free_poly(&p); - free_fft_settings(&fs); - free_kzg_settings(&ks); - free_fk20_single_settings(&fk); -} - -void fk_multi_settings(void) { - FFTSettings fs; - KZGSettings ks; - FK20MultiSettings fk; - uint64_t n = 5; - uint64_t secrets_len = 33; - g1_t s1[secrets_len]; - g2_t s2[secrets_len]; - - // Initialise the secrets and data structures - generate_trusted_setup(s1, s2, &secret, secrets_len); - TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, n)); - TEST_CHECK(C_KZG_OK == new_kzg_settings(&ks, s1, s2, secrets_len, &fs)); - TEST_CHECK(C_KZG_OK == new_fk20_multi_settings(&fk, 32, 4, &ks)); - - // Don't do anything. Run this with `valgrind` to check that memory is correctly allocated and freed. - - free_fft_settings(&fs); - free_kzg_settings(&ks); - free_fk20_multi_settings(&fk); -} - -void fk_multi_0(void) { - FFTSettings fs; - KZGSettings ks; - FK20MultiSettings fk; - uint64_t n, chunk_len, chunk_count; - uint64_t secrets_len; - g1_t *s1; - g2_t *s2; - poly p; - uint64_t vv[] = {1, 2, 3, 4, 7, 8, 9, 10, 13, 14, 1, 15, 1, 1000, 134, 33}; - g1_t commitment; - g1_t *all_proofs; - fr_t *extended_coeffs, *extended_coeffs_fft; - fr_t *ys, *ys2; - uint64_t domain_stride; - - chunk_len = 16; - chunk_count = 32; - n = chunk_len * chunk_count; - secrets_len = 2 * n; - - TEST_CHECK(C_KZG_OK == new_g1_array(&s1, secrets_len)); - TEST_CHECK(C_KZG_OK == new_g2_array(&s2, secrets_len)); - - generate_trusted_setup(s1, s2, &secret, secrets_len); - TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 4 + 5 + 1)); - TEST_CHECK(C_KZG_OK == new_kzg_settings(&ks, s1, s2, secrets_len, &fs)); - TEST_CHECK(C_KZG_OK == new_fk20_multi_settings(&fk, n * 2, chunk_len, &ks)); - - // Create a test polynomial: 512 coefficients - TEST_CHECK(C_KZG_OK == new_poly(&p, n)); - for (int i = 0; i < chunk_count; i++) { - for (int j = 0; j < chunk_len; j++) { - uint64_t v = vv[j]; - if (j == 3) v += i; - if (j == 5) v += i * i; - fr_from_uint64(&p.coeffs[i * chunk_len + j], v); - } - fr_negate(&p.coeffs[i * chunk_len + 12], &p.coeffs[i * chunk_len + 12]); - fr_negate(&p.coeffs[i * chunk_len + 14], &p.coeffs[i * chunk_len + 14]); - } - - TEST_CHECK(C_KZG_OK == commit_to_poly(&commitment, &p, &ks)); - - // Compute the multi proofs, assuming that the polynomial will be extended with zeros - TEST_CHECK(C_KZG_OK == new_g1_array(&all_proofs, 2 * chunk_count)); - TEST_CHECK(C_KZG_OK == da_using_fk20_multi(all_proofs, &p, &fk)); - - // Now actually extend the polynomial with zeros - TEST_CHECK(C_KZG_OK == new_fr_array(&extended_coeffs, 2 * n)); - for (uint64_t i = 0; i < n; i++) { - extended_coeffs[i] = p.coeffs[i]; - } - for (uint64_t i = n; i < 2 * n; i++) { - extended_coeffs[i] = fr_zero; - } - TEST_CHECK(C_KZG_OK == new_fr_array(&extended_coeffs_fft, 2 * n)); - TEST_CHECK(C_KZG_OK == fft_fr(extended_coeffs_fft, extended_coeffs, false, 2 * n, &fs)); - TEST_CHECK(C_KZG_OK == reverse_bit_order(extended_coeffs_fft, sizeof extended_coeffs_fft[0], 2 * n)); - - // Verify the proofs - TEST_CHECK(C_KZG_OK == new_fr_array(&ys, chunk_len)); - TEST_CHECK(C_KZG_OK == new_fr_array(&ys2, chunk_len)); - domain_stride = fs.max_width / (2 * n); - for (uint64_t pos = 0; pos < 2 * chunk_count; pos++) { - uint64_t domain_pos, stride; - fr_t x; - bool result; - - domain_pos = reverse_bits_limited(2 * chunk_count, pos); - x = fs.expanded_roots_of_unity[domain_pos * domain_stride]; - - // The ys from the extended coeffients - for (uint64_t i = 0; i < chunk_len; i++) { - ys[i] = extended_coeffs_fft[chunk_len * pos + i]; - } - TEST_CHECK(C_KZG_OK == reverse_bit_order(ys, sizeof ys[0], chunk_len)); - - // Now recreate the ys by evaluating the polynomial in the sub-domain range - stride = fs.max_width / chunk_len; - for (uint64_t i = 0; i < chunk_len; i++) { - fr_t z; - fr_mul(&z, &x, &fs.expanded_roots_of_unity[i * stride]); - eval_poly(&ys2[i], &p, &z); - } - - // ys and ys2 should be equal - for (uint64_t i = 0; i < chunk_len; i++) { - TEST_CHECK(fr_equal(&ys[i], &ys2[i])); - } - - // Verify this proof - TEST_CHECK(C_KZG_OK == check_proof_multi(&result, &commitment, &all_proofs[pos], &x, ys, chunk_len, &ks)); - TEST_CHECK(true == result); - } - - free_poly(&p); - free(all_proofs); - free(extended_coeffs); - free(extended_coeffs_fft); - free(ys); - free(ys2); - free(s1); - free(s2); - free_fft_settings(&fs); - free_kzg_settings(&ks); - free_fk20_multi_settings(&fk); -} - -// TODO: compare results of fk20_multi_da_opt() and fk20_compute_proof_multi() - -TEST_LIST = { - {"FK20_PROOFS_TEST", title}, - {"fk_single", fk_single}, - {"fk_single_strided", fk_single_strided}, - {"fk_multi_settings", fk_multi_settings}, - {"fk_multi_0", fk_multi_0}, - {NULL, NULL} /* zero record marks the end of the list */ -}; diff --git a/src/kzg_proofs.c b/src/kzg_proofs.c index 7216e0b..e640d7e 100644 --- a/src/kzg_proofs.c +++ b/src/kzg_proofs.c @@ -239,3 +239,172 @@ void free_kzg_settings(KZGSettings *ks) { free(ks->secret_g2); ks->length = 0; } + +#ifdef KZGTEST + +#include "../inc/acutest.h" +#include "test_util.h" + +void proof_single(void) { + // Our polynomial: degree 15, 16 coefficients + uint64_t coeffs[] = {1, 2, 3, 4, 7, 7, 7, 7, 13, 13, 13, 13, 13, 13, 13, 13}; + int poly_len = sizeof coeffs / sizeof coeffs[0]; + uint64_t secrets_len = poly_len + 1; + + FFTSettings fs; + KZGSettings ks; + g1_t s1[secrets_len]; + g2_t s2[secrets_len]; + poly p; + g1_t commitment, proof; + fr_t x, value; + bool result; + + // Create the polynomial + new_poly(&p, poly_len); + for (int i = 0; i < poly_len; i++) { + fr_from_uint64(&p.coeffs[i], coeffs[i]); + } + + // Initialise the secrets and data structures + generate_trusted_setup(s1, s2, &secret, secrets_len); + TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 4)); // ln_2 of poly_len + TEST_CHECK(C_KZG_OK == new_kzg_settings(&ks, s1, s2, secrets_len, &fs)); + + // Compute the proof for x = 25 + fr_from_uint64(&x, 25); + TEST_CHECK(C_KZG_OK == commit_to_poly(&commitment, &p, &ks)); + TEST_CHECK(C_KZG_OK == compute_proof_single(&proof, &p, &x, &ks)); + + eval_poly(&value, &p, &x); + + // Verify the proof that the (unknown) polynomial has y = value at x = 25 + TEST_CHECK(C_KZG_OK == check_proof_single(&result, &commitment, &proof, &x, &value, &ks)); + TEST_CHECK(true == result); + + // Change the value and check that the proof fails + fr_add(&value, &value, &fr_one); + TEST_CHECK(C_KZG_OK == check_proof_single(&result, &commitment, &proof, &x, &value, &ks)); + TEST_CHECK(false == result); + + free_fft_settings(&fs); + free_kzg_settings(&ks); + free_poly(&p); +} + +void proof_multi(void) { + // Our polynomial: degree 15, 16 coefficients + uint64_t coeffs[] = {1, 2, 3, 4, 7, 7, 7, 7, 13, 13, 13, 13, 13, 13, 13, 13}; + int poly_len = sizeof coeffs / sizeof coeffs[0]; + + FFTSettings fs1, fs2; + KZGSettings ks1, ks2; + poly p; + g1_t commitment, proof; + fr_t x, tmp; + bool result; + + // Compute proof at 2^coset_scale points + int coset_scale = 3, coset_len = (1 << coset_scale); + fr_t y[coset_len]; + + uint64_t secrets_len = poly_len > coset_len ? poly_len + 1 : coset_len + 1; + g1_t s1[secrets_len]; + g2_t s2[secrets_len]; + + // Create the polynomial + new_poly(&p, poly_len); + for (int i = 0; i < poly_len; i++) { + fr_from_uint64(&p.coeffs[i], coeffs[i]); + } + + // Initialise the secrets and data structures + generate_trusted_setup(s1, s2, &secret, secrets_len); + TEST_CHECK(C_KZG_OK == new_fft_settings(&fs1, 4)); // ln_2 of poly_len + TEST_CHECK(C_KZG_OK == new_kzg_settings(&ks1, s1, s2, secrets_len, &fs1)); + + // Commit to the polynomial + TEST_CHECK(C_KZG_OK == commit_to_poly(&commitment, &p, &ks1)); + + TEST_CHECK(C_KZG_OK == new_fft_settings(&fs2, coset_scale)); + TEST_CHECK(C_KZG_OK == new_kzg_settings(&ks2, s1, s2, secrets_len, &fs2)); + + // Compute proof at the points [x * root_i] 0 <= i < coset_len + fr_from_uint64(&x, 5431); + TEST_CHECK(C_KZG_OK == compute_proof_multi(&proof, &p, &x, coset_len, &ks2)); + + // y_i is the value of the polynomial at each x_i + for (int i = 0; i < coset_len; i++) { + fr_mul(&tmp, &x, &ks2.fs->expanded_roots_of_unity[i]); + eval_poly(&y[i], &p, &tmp); + } + + // Verify the proof that the (unknown) polynomial has value y_i at x_i + TEST_CHECK(C_KZG_OK == check_proof_multi(&result, &commitment, &proof, &x, y, coset_len, &ks2)); + TEST_CHECK(true == result); + + // Change a value and check that the proof fails + fr_add(y + coset_len / 2, y + coset_len / 2, &fr_one); + TEST_CHECK(C_KZG_OK == check_proof_multi(&result, &commitment, &proof, &x, y, coset_len, &ks2)); + TEST_CHECK(false == result); + + free_fft_settings(&fs1); + free_fft_settings(&fs2); + free_kzg_settings(&ks1); + free_kzg_settings(&ks2); + free_poly(&p); +} + +void commit_to_nil_poly(void) { + poly a; + FFTSettings fs; + KZGSettings ks; + uint64_t secrets_len = 16; + g1_t s1[secrets_len]; + g2_t s2[secrets_len]; + g1_t result; + + // Initialise the (arbitrary) secrets and data structures + generate_trusted_setup(s1, s2, &secret, secrets_len); + TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 4)); + TEST_CHECK(C_KZG_OK == new_kzg_settings(&ks, s1, s2, secrets_len, &fs)); + + new_poly(&a, 0); + TEST_CHECK(C_KZG_OK == commit_to_poly(&result, &a, &ks)); + TEST_CHECK(g1_equal(&g1_identity, &result)); + + free_fft_settings(&fs); + free_kzg_settings(&ks); +} + +void commit_to_too_long_poly(void) { + poly a; + FFTSettings fs; + KZGSettings ks; + uint64_t poly_len = 32, secrets_len = 16; // poly is longer than secrets! + g1_t s1[secrets_len]; + g2_t s2[secrets_len]; + g1_t result; + + // Initialise the (arbitrary) secrets and data structures + generate_trusted_setup(s1, s2, &secret, secrets_len); + TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 4)); + TEST_CHECK(C_KZG_OK == new_kzg_settings(&ks, s1, s2, secrets_len, &fs)); + + new_poly(&a, poly_len); + TEST_CHECK(C_KZG_BADARGS == commit_to_poly(&result, &a, &ks)); + + free_fft_settings(&fs); + free_kzg_settings(&ks); +} + +TEST_LIST = { + {"KZG_PROOFS_TEST", title}, + {"proof_single", proof_single}, + {"proof_multi", proof_multi}, + {"commit_to_nil_poly", commit_to_nil_poly}, + {"commit_to_too_long_poly", commit_to_too_long_poly}, + {NULL, NULL} /* zero record marks the end of the list */ +}; + +#endif // KZGTEST diff --git a/src/kzg_proofs_test.c b/src/kzg_proofs_test.c deleted file mode 100644 index 83c0267..0000000 --- a/src/kzg_proofs_test.c +++ /dev/null @@ -1,181 +0,0 @@ -/* - * Copyright 2021 Benjamin Edgington - * - * Licensed under the Apache License, Version 2.0 (the "License"); - * you may not use this file except in compliance with the License. - * You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -#include "../inc/acutest.h" -#include "test_util.h" -#include "kzg_proofs.h" - -void proof_single(void) { - // Our polynomial: degree 15, 16 coefficients - uint64_t coeffs[] = {1, 2, 3, 4, 7, 7, 7, 7, 13, 13, 13, 13, 13, 13, 13, 13}; - int poly_len = sizeof coeffs / sizeof coeffs[0]; - uint64_t secrets_len = poly_len + 1; - - FFTSettings fs; - KZGSettings ks; - g1_t s1[secrets_len]; - g2_t s2[secrets_len]; - poly p; - g1_t commitment, proof; - fr_t x, value; - bool result; - - // Create the polynomial - new_poly(&p, poly_len); - for (int i = 0; i < poly_len; i++) { - fr_from_uint64(&p.coeffs[i], coeffs[i]); - } - - // Initialise the secrets and data structures - generate_trusted_setup(s1, s2, &secret, secrets_len); - TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 4)); // ln_2 of poly_len - TEST_CHECK(C_KZG_OK == new_kzg_settings(&ks, s1, s2, secrets_len, &fs)); - - // Compute the proof for x = 25 - fr_from_uint64(&x, 25); - TEST_CHECK(C_KZG_OK == commit_to_poly(&commitment, &p, &ks)); - TEST_CHECK(C_KZG_OK == compute_proof_single(&proof, &p, &x, &ks)); - - eval_poly(&value, &p, &x); - - // Verify the proof that the (unknown) polynomial has y = value at x = 25 - TEST_CHECK(C_KZG_OK == check_proof_single(&result, &commitment, &proof, &x, &value, &ks)); - TEST_CHECK(true == result); - - // Change the value and check that the proof fails - fr_add(&value, &value, &fr_one); - TEST_CHECK(C_KZG_OK == check_proof_single(&result, &commitment, &proof, &x, &value, &ks)); - TEST_CHECK(false == result); - - free_fft_settings(&fs); - free_kzg_settings(&ks); - free_poly(&p); -} - -void proof_multi(void) { - // Our polynomial: degree 15, 16 coefficients - uint64_t coeffs[] = {1, 2, 3, 4, 7, 7, 7, 7, 13, 13, 13, 13, 13, 13, 13, 13}; - int poly_len = sizeof coeffs / sizeof coeffs[0]; - - FFTSettings fs1, fs2; - KZGSettings ks1, ks2; - poly p; - g1_t commitment, proof; - fr_t x, tmp; - bool result; - - // Compute proof at 2^coset_scale points - int coset_scale = 3, coset_len = (1 << coset_scale); - fr_t y[coset_len]; - - uint64_t secrets_len = poly_len > coset_len ? poly_len + 1 : coset_len + 1; - g1_t s1[secrets_len]; - g2_t s2[secrets_len]; - - // Create the polynomial - new_poly(&p, poly_len); - for (int i = 0; i < poly_len; i++) { - fr_from_uint64(&p.coeffs[i], coeffs[i]); - } - - // Initialise the secrets and data structures - generate_trusted_setup(s1, s2, &secret, secrets_len); - TEST_CHECK(C_KZG_OK == new_fft_settings(&fs1, 4)); // ln_2 of poly_len - TEST_CHECK(C_KZG_OK == new_kzg_settings(&ks1, s1, s2, secrets_len, &fs1)); - - // Commit to the polynomial - TEST_CHECK(C_KZG_OK == commit_to_poly(&commitment, &p, &ks1)); - - TEST_CHECK(C_KZG_OK == new_fft_settings(&fs2, coset_scale)); - TEST_CHECK(C_KZG_OK == new_kzg_settings(&ks2, s1, s2, secrets_len, &fs2)); - - // Compute proof at the points [x * root_i] 0 <= i < coset_len - fr_from_uint64(&x, 5431); - TEST_CHECK(C_KZG_OK == compute_proof_multi(&proof, &p, &x, coset_len, &ks2)); - - // y_i is the value of the polynomial at each x_i - for (int i = 0; i < coset_len; i++) { - fr_mul(&tmp, &x, &ks2.fs->expanded_roots_of_unity[i]); - eval_poly(&y[i], &p, &tmp); - } - - // Verify the proof that the (unknown) polynomial has value y_i at x_i - TEST_CHECK(C_KZG_OK == check_proof_multi(&result, &commitment, &proof, &x, y, coset_len, &ks2)); - TEST_CHECK(true == result); - - // Change a value and check that the proof fails - fr_add(y + coset_len / 2, y + coset_len / 2, &fr_one); - TEST_CHECK(C_KZG_OK == check_proof_multi(&result, &commitment, &proof, &x, y, coset_len, &ks2)); - TEST_CHECK(false == result); - - free_fft_settings(&fs1); - free_fft_settings(&fs2); - free_kzg_settings(&ks1); - free_kzg_settings(&ks2); - free_poly(&p); -} - -void commit_to_nil_poly(void) { - poly a; - FFTSettings fs; - KZGSettings ks; - uint64_t secrets_len = 16; - g1_t s1[secrets_len]; - g2_t s2[secrets_len]; - g1_t result; - - // Initialise the (arbitrary) secrets and data structures - generate_trusted_setup(s1, s2, &secret, secrets_len); - TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 4)); - TEST_CHECK(C_KZG_OK == new_kzg_settings(&ks, s1, s2, secrets_len, &fs)); - - new_poly(&a, 0); - TEST_CHECK(C_KZG_OK == commit_to_poly(&result, &a, &ks)); - TEST_CHECK(g1_equal(&g1_identity, &result)); - - free_fft_settings(&fs); - free_kzg_settings(&ks); -} - -void commit_to_too_long_poly(void) { - poly a; - FFTSettings fs; - KZGSettings ks; - uint64_t poly_len = 32, secrets_len = 16; // poly is longer than secrets! - g1_t s1[secrets_len]; - g2_t s2[secrets_len]; - g1_t result; - - // Initialise the (arbitrary) secrets and data structures - generate_trusted_setup(s1, s2, &secret, secrets_len); - TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 4)); - TEST_CHECK(C_KZG_OK == new_kzg_settings(&ks, s1, s2, secrets_len, &fs)); - - new_poly(&a, poly_len); - TEST_CHECK(C_KZG_BADARGS == commit_to_poly(&result, &a, &ks)); - - free_fft_settings(&fs); - free_kzg_settings(&ks); -} - -TEST_LIST = { - {"KZG_PROOFS_TEST", title}, - {"proof_single", proof_single}, - {"proof_multi", proof_multi}, - {"commit_to_nil_poly", commit_to_nil_poly}, - {"commit_to_too_long_poly", commit_to_too_long_poly}, - {NULL, NULL} /* zero record marks the end of the list */ -}; diff --git a/src/poly.c b/src/poly.c index 73d0dc4..fee2ba5 100644 --- a/src/poly.c +++ b/src/poly.c @@ -35,6 +35,24 @@ static uint64_t poly_quotient_length(const poly *dividend, const poly *divisor) return dividend->length >= divisor->length ? dividend->length - divisor->length + 1 : 0; } +/** + * Pad with zeros or truncate an array of field elements to a specific size. + * + * @param[out] out The padded/truncated array + * @param[in] in The original array to be padded + * @param[in] n_in The number of elements of @p in to take + * @param[in] n_out The length of @p out + */ +static void pad(fr_t *out, const fr_t *in, uint64_t n_in, uint64_t n_out) { + uint64_t num = min_u64(n_in, n_out); + for (int i = 0; i < num; i++) { + out[i] = in[i]; + } + for (int i = num; i < n_out; i++) { + out[i] = fr_zero; + } +} + /** * Return a copy of a polynomial ensuring that the order is correct. * @@ -102,7 +120,7 @@ void eval_poly(fr_t *out, const poly *p, const fr_t *x) { * @retval C_CZK_BADARGS Invalid parameters were supplied * @retval C_CZK_MALLOC Memory allocation failed */ -C_KZG_RET poly_long_div(poly *out, const poly *dividend, const poly *divisor) { +static C_KZG_RET poly_long_div(poly *out, const poly *dividend, const poly *divisor) { uint64_t a_pos = dividend->length - 1; uint64_t b_pos = divisor->length - 1; uint64_t diff = a_pos - b_pos; @@ -156,7 +174,7 @@ C_KZG_RET poly_long_div(poly *out, const poly *dividend, const poly *divisor) { * @param[in] b The multiplier polynomial * @retval C_CZK_OK All is well */ -C_KZG_RET poly_mul_direct(poly *out, const poly *a, const poly *b) { +static C_KZG_RET poly_mul_direct(poly *out, const poly *a, const poly *b) { uint64_t a_degree = a->length - 1; uint64_t b_degree = b->length - 1; @@ -176,23 +194,6 @@ C_KZG_RET poly_mul_direct(poly *out, const poly *a, const poly *b) { return C_KZG_OK; } -/** - * Pad with zeros or truncate an array of field elements to a specific size. - * - * @param[out] out The padded/truncated array - * @param[in] in The original array to be padded - * @param[in] n_in The number of elements of @p in to take - * @param[in] n_out The length of @p out - */ -void pad(fr_t *out, const fr_t *in, uint64_t n_in, uint64_t n_out) { - uint64_t num = min_u64(n_in, n_out); - for (int i = 0; i < num; i++) { - out[i] = in[i]; - } - for (int i = num; i < n_out; i++) { - out[i] = fr_zero; - } -} /** * Calculate the (possibly truncated) product of two polynomials. @@ -212,7 +213,7 @@ void pad(fr_t *out, const fr_t *in, uint64_t n_in, uint64_t n_out) { * @retval C_CZK_ERROR An internal error occurred * @retval C_CZK_MALLOC Memory allocation failed */ -C_KZG_RET poly_mul_fft(poly *out, const poly *a, const poly *b, FFTSettings *fs_) { +static C_KZG_RET poly_mul_fft(poly *out, const poly *a, const poly *b, FFTSettings *fs_) { // Truncate a and b so as not to do excess work for the number of coefficients required. uint64_t a_len = min_u64(a->length, out->length); @@ -395,7 +396,7 @@ C_KZG_RET poly_flip(poly *out, const poly *in) { * @retval C_CZK_ERROR An internal error occurred * @retval C_CZK_MALLOC Memory allocation failed */ -C_KZG_RET poly_fast_div(poly *out, const poly *dividend, const poly *divisor) { +static C_KZG_RET poly_fast_div(poly *out, const poly *dividend, const poly *divisor) { // Dividing by zero is undefined CHECK(divisor->length > 0); @@ -467,8 +468,6 @@ C_KZG_RET poly_fast_div(poly *out, const poly *dividend, const poly *divisor) { * @retval C_CZK_BADARGS Invalid parameters were supplied * @retval C_CZK_ERROR An internal error occurred * @retval C_CZK_MALLOC Memory allocation failed - * - * @todo Implement routine selection */ C_KZG_RET poly_mul_(poly *out, const poly *a, const poly *b, FFTSettings *fs) { if (a->length < 64 || b->length < 64 || out->length < 128) { // Tunable parameter @@ -491,8 +490,6 @@ C_KZG_RET poly_mul_(poly *out, const poly *a, const poly *b, FFTSettings *fs) { * @retval C_CZK_BADARGS Invalid parameters were supplied * @retval C_CZK_ERROR An internal error occurred * @retval C_CZK_MALLOC Memory allocation failed - * - * @todo Implement routine selection */ C_KZG_RET poly_mul(poly *out, const poly *a, const poly *b) { return poly_mul_(out, a, b, NULL); @@ -578,3 +575,407 @@ void free_poly(poly *p) { free(p->coeffs); } } + +#ifdef KZGTEST + +#include "../inc/acutest.h" +#include "test_util.h" + +// If anyone has a nicer way to initialise these test data, I'd love to hear it. +typedef struct polydata { + int length; + int coeffs[]; +} polydata; + +// (x^2 - 1) / (x + 1) = x - 1 +polydata test_0_0 = {3, {-1, 0, 1}}; +polydata test_0_1 = {2, {1, 1}}; +polydata test_0_2 = {2, {-1, 1}}; + +// (12x^3 - 11x^2 + 9x + 18) / (4x + 3) = 3x^2 - 5x + 6 +polydata test_1_0 = {4, {18, 9, -11, 12}}; +polydata test_1_1 = {2, {3, 4}}; +polydata test_1_2 = {3, {6, -5, 3}}; + +// (x + 1) / (x^2 - 1) = nil +polydata test_2_0 = {2, {1, 1}}; +polydata test_2_1 = {3, {-1, 0, 2}}; +polydata test_2_2 = {0, {}}; + +// (10x^2 + 20x + 30) / 10 = x^2 + 2x + 3 +polydata test_3_0 = {3, {30, 20, 10}}; +polydata test_3_1 = {1, {10}}; +polydata test_3_2 = {3, {3, 2, 1}}; + +// (x^2 + x) / (x + 1) = x +polydata test_4_0 = {3, {0, 1, 1}}; +polydata test_4_1 = {2, {1, 1}}; +polydata test_4_2 = {2, {0, 1}}; + +// (x^2 + x + 1) / 1 = x^2 + x + 1 +polydata test_5_0 = {3, {1, 1, 1}}; +polydata test_5_1 = {1, {1}}; +polydata test_5_2 = {3, {1, 1, 1}}; + +// (x^2 + x + 1) / (0x + 1) = x^2 + x + 1 +polydata test_6_0 = {3, {1, 1, 1}}; +polydata test_6_1 = {2, {1, 0}}; // The highest coefficient is zero +polydata test_6_2 = {3, {1, 1, 1}}; + +polydata *test[][3] = {{&test_0_0, &test_0_1, &test_0_2}, {&test_1_0, &test_1_1, &test_1_2}, + {&test_2_0, &test_2_1, &test_2_2}, {&test_3_0, &test_3_1, &test_3_2}, + {&test_4_0, &test_4_1, &test_4_2}, {&test_5_0, &test_5_1, &test_5_2}, + {&test_6_0, &test_6_1, &test_6_2}}; + +/* Internal utility function */ +void new_test_poly(poly *p, polydata *data) { + new_poly(p, data->length); + for (int i = 0; i < p->length; i++) { + int coeff = data->coeffs[i]; + if (coeff >= 0) { + fr_from_uint64(&p->coeffs[i], coeff); + } else { + fr_from_uint64(&p->coeffs[i], -coeff); + fr_negate(&p->coeffs[i], &p->coeffs[i]); + } + } +} + +void poly_test_div(void) { + poly dividend, divisor, expected, actual; + int ntest = sizeof test / sizeof test[0]; + + for (int i = 0; i < ntest; i++) { + new_test_poly(÷nd, test[i][0]); + new_test_poly(&divisor, test[i][1]); + new_test_poly(&expected, test[i][2]); + + if (TEST_CHECK(C_KZG_OK == new_poly_div(&actual, ÷nd, &divisor))) { + if (TEST_CHECK(actual.length == expected.length)) { + for (int j = 0; j < actual.length; j++) { + TEST_CHECK(fr_equal(&actual.coeffs[j], &expected.coeffs[j])); + TEST_MSG("Failed test %d with incorrect value", i); + } + } else { + TEST_MSG("Failed test %d with incorrect length.", i); + } + } else { + TEST_MSG("Failed test %d with bad return value.", i); + } + + free_poly(÷nd); + free_poly(&divisor); + free_poly(&expected); + free_poly(&actual); + } +} + +void poly_div_by_zero(void) { + fr_t a[2]; + poly dividend, divisor, dummy; + + // Calculate (x + 1) / 0 = FAIL + + // Dividend + fr_from_uint64(&a[0], 1); + fr_from_uint64(&a[1], 1); + dividend.length = 2; + dividend.coeffs = a; + + // Divisor + new_poly(&divisor, 0); + + TEST_CHECK(C_KZG_BADARGS == new_poly_div(&dummy, ÷nd, &divisor)); + + free_poly(&divisor); + free_poly(&dummy); +} + +void poly_eval_check(void) { + uint64_t n = 10; + fr_t actual, expected; + poly p; + new_poly(&p, n); + for (uint64_t i = 0; i < n; i++) { + fr_from_uint64(&p.coeffs[i], i + 1); + } + fr_from_uint64(&expected, n * (n + 1) / 2); + + eval_poly(&actual, &p, &fr_one); + + TEST_CHECK(fr_equal(&expected, &actual)); + + free_poly(&p); +} + +void poly_eval_0_check(void) { + uint64_t n = 7, a = 597; + fr_t actual, expected; + poly p; + new_poly(&p, n); + for (uint64_t i = 0; i < n; i++) { + fr_from_uint64(&p.coeffs[i], i + a); + } + fr_from_uint64(&expected, a); + + eval_poly(&actual, &p, &fr_zero); + + TEST_CHECK(fr_equal(&expected, &actual)); + + free_poly(&p); +} + +void poly_eval_nil_check(void) { + uint64_t n = 0; + fr_t actual; + poly p; + new_poly(&p, n); + + eval_poly(&actual, &p, &fr_one); + + TEST_CHECK(fr_equal(&fr_zero, &actual)); + + free_poly(&p); +} + +void poly_mul_direct_test(void) { + + // Calculate (3x^2 - 5x + 6) * (4x + 3) = 12x^3 - 11x^2 + 9x + 18 + static polydata multiplier_data = {3, {6, -5, 3}}; + static polydata multiplicand_data = {2, {3, 4}}; + static polydata expected_data = {4, {18, 9, -11, 12}}; + + poly multiplicand, multiplier, expected, actual0, actual1; + + new_test_poly(&multiplicand, &multiplicand_data); + new_test_poly(&multiplier, &multiplier_data); + new_test_poly(&expected, &expected_data); + + new_poly(&actual0, 4); + TEST_CHECK(C_KZG_OK == poly_mul_direct(&actual0, &multiplicand, &multiplier)); + TEST_CHECK(fr_equal(&expected.coeffs[0], &actual0.coeffs[0])); + TEST_CHECK(fr_equal(&expected.coeffs[1], &actual0.coeffs[1])); + TEST_CHECK(fr_equal(&expected.coeffs[2], &actual0.coeffs[2])); + TEST_CHECK(fr_equal(&expected.coeffs[3], &actual0.coeffs[3])); + + // Check commutativity + new_poly(&actual1, 4); + TEST_CHECK(C_KZG_OK == poly_mul_direct(&actual1, &multiplier, &multiplicand)); + TEST_CHECK(fr_equal(&expected.coeffs[0], &actual1.coeffs[0])); + TEST_CHECK(fr_equal(&expected.coeffs[1], &actual1.coeffs[1])); + TEST_CHECK(fr_equal(&expected.coeffs[2], &actual1.coeffs[2])); + TEST_CHECK(fr_equal(&expected.coeffs[3], &actual1.coeffs[3])); + + free_poly(&multiplicand); + free_poly(&multiplier); + free_poly(&expected); + free_poly(&actual0); + free_poly(&actual1); +} + +void poly_mul_fft_test(void) { + + // Calculate (3x^2 - 5x + 6) * (4x + 3) = 12x^3 - 11x^2 + 9x + 18 + static polydata multiplier_data = {3, {6, -5, 3}}; + static polydata multiplicand_data = {2, {3, 4}}; + static polydata expected_data = {4, {18, 9, -11, 12}}; + + poly multiplicand, multiplier, expected, actual0, actual1; + + new_test_poly(&multiplicand, &multiplicand_data); + new_test_poly(&multiplier, &multiplier_data); + new_test_poly(&expected, &expected_data); + + new_poly(&actual0, 4); + TEST_CHECK(C_KZG_OK == poly_mul_fft(&actual0, &multiplicand, &multiplier, NULL)); + TEST_CHECK(fr_equal(&expected.coeffs[0], &actual0.coeffs[0])); + TEST_CHECK(fr_equal(&expected.coeffs[1], &actual0.coeffs[1])); + TEST_CHECK(fr_equal(&expected.coeffs[2], &actual0.coeffs[2])); + TEST_CHECK(fr_equal(&expected.coeffs[3], &actual0.coeffs[3])); + + // Check commutativity + new_poly(&actual1, 4); + TEST_CHECK(C_KZG_OK == poly_mul_fft(&actual1, &multiplier, &multiplicand, NULL)); + TEST_CHECK(fr_equal(&expected.coeffs[0], &actual1.coeffs[0])); + TEST_CHECK(fr_equal(&expected.coeffs[1], &actual1.coeffs[1])); + TEST_CHECK(fr_equal(&expected.coeffs[2], &actual1.coeffs[2])); + TEST_CHECK(fr_equal(&expected.coeffs[3], &actual1.coeffs[3])); + + free_poly(&multiplicand); + free_poly(&multiplier); + free_poly(&expected); + free_poly(&actual0); + free_poly(&actual1); +} + +void poly_inverse_simple_0(void) { + + // 1 / (1 - x) = 1 + x + x^2 + ... + + poly p, q; + int d = 16; // The number of terms to take + + new_poly(&p, 2); + p.coeffs[0] = fr_one; + p.coeffs[1] = fr_one; + fr_negate(&p.coeffs[1], &p.coeffs[1]); + + new_poly(&q, d); + TEST_CHECK(C_KZG_OK == poly_inverse(&q, &p)); + + for (int i = 0; i < d; i++) { + TEST_CHECK(fr_is_one(&q.coeffs[i])); + } + + free_poly(&p); + free_poly(&q); +} + +void poly_inverse_simple_1(void) { + + // 1 / (1 + x) = 1 - x + x^2 - ... + + poly p, q; + int d = 16; // The number of terms to take + + new_poly(&p, 2); + p.coeffs[0] = fr_one; + p.coeffs[1] = fr_one; + + new_poly(&q, d); + TEST_CHECK(C_KZG_OK == poly_inverse(&q, &p)); + + for (int i = 0; i < d; i++) { + fr_t tmp = q.coeffs[i]; + if (i & 1) { + fr_negate(&tmp, &tmp); + } + TEST_CHECK(fr_is_one(&tmp)); + } + + free_poly(&p); + free_poly(&q); +} + +void poly_mul_random(void) { + + // Compare the output of poly_mul_direct() and poly_mul_fft() + + poly multiplicand, multiplier; + poly q0, q1; + + for (int k = 0; k < 256; k++) { + + int multiplicand_length = 1 + rand() % 1000; + int multiplier_length = 1 + rand() % 1000; + int out_length = 1 + rand() % 1000; + + new_poly(&multiplicand, multiplicand_length); + new_poly(&multiplier, multiplier_length); + + for (int i = 0; i < multiplicand_length; i++) { + multiplicand.coeffs[i] = rand_fr(); + } + for (int i = 0; i < multiplier_length; i++) { + multiplier.coeffs[i] = rand_fr(); + } + + // Ensure that the polynomials' orders corresponds to their lengths + if (fr_is_zero(&multiplicand.coeffs[multiplicand.length - 1])) { + multiplicand.coeffs[multiplicand.length - 1] = fr_one; + } + if (fr_is_zero(&multiplier.coeffs[multiplier.length - 1])) { + multiplier.coeffs[multiplier.length - 1] = fr_one; + } + + new_poly(&q0, out_length); // Truncate the result + TEST_CHECK(C_KZG_OK == poly_mul_direct(&q0, &multiplicand, &multiplier)); + + new_poly(&q1, out_length); + TEST_CHECK(C_KZG_OK == poly_mul_fft(&q1, &multiplicand, &multiplier, NULL)); + + TEST_CHECK(q1.length == q0.length); + for (int i = 0; i < q0.length; i++) { + if (!TEST_CHECK(fr_equal(&q0.coeffs[i], &q1.coeffs[i]))) { + TEST_MSG("round = %d, i = %d, multiplicand_length = %lu, multiplier_length = %lu, out_length = %lu", k, + i, multiplicand.length, multiplier.length, q0.length); + } + } + + free_poly(&multiplicand); + free_poly(&multiplier); + free_poly(&q0); + free_poly(&q1); + } +} + +void poly_div_random(void) { + + // Compare the output of poly_fast_div() and poly_long_div() + + poly dividend, divisor; + poly q0, q1; + + for (int k = 0; k < 256; k++) { + + int dividend_length = 2 + rand() % 1000; + int divisor_length = 1 + rand() % dividend_length; + + new_poly(÷nd, dividend_length); + new_poly(&divisor, divisor_length); + + for (int i = 0; i < dividend_length; i++) { + dividend.coeffs[i] = rand_fr(); + } + for (int i = 0; i < divisor_length; i++) { + divisor.coeffs[i] = rand_fr(); + } + + // Ensure that the polynomials' orders corresponds to their lengths + if (fr_is_zero(÷nd.coeffs[dividend.length - 1])) { + dividend.coeffs[dividend.length - 1] = fr_one; + } + if (fr_is_zero(&divisor.coeffs[divisor.length - 1])) { + divisor.coeffs[divisor.length - 1] = fr_one; + } + + new_poly(&q0, dividend.length - divisor.length + 1); + TEST_CHECK(C_KZG_OK == poly_long_div(&q0, ÷nd, &divisor)); + + new_poly(&q1, dividend.length - divisor.length + 1); + TEST_CHECK(C_KZG_OK == poly_fast_div(&q1, ÷nd, &divisor)); + + TEST_CHECK(q1.length == q0.length); + for (int i = 0; i < q0.length; i++) { + if (!TEST_CHECK(fr_equal(&q0.coeffs[i], &q1.coeffs[i]))) { + TEST_MSG("round = %d, dividend_length = %lu, divisor_length = %lu, i = %d", k, dividend.length, + divisor.length, i); + } + } + + free_poly(÷nd); + free_poly(&divisor); + free_poly(&q0); + free_poly(&q1); + } +} + +TEST_LIST = { + {"POLY_TEST", title}, + {"poly_test_div", poly_test_div}, +#ifndef DEBUG + {"poly_div_by_zero", poly_div_by_zero}, +#endif + {"poly_eval_check", poly_eval_check}, + {"poly_eval_0_check", poly_eval_0_check}, + {"poly_eval_nil_check", poly_eval_nil_check}, + {"poly_mul_direct_test", poly_mul_direct_test}, + {"poly_mul_fft_test", poly_mul_fft_test}, + {"poly_inverse_simple_0", poly_inverse_simple_0}, + {"poly_inverse_simple_1", poly_inverse_simple_1}, + {"poly_mul_random", poly_mul_random}, + {"poly_div_random", poly_div_random}, + {NULL, NULL} /* zero record marks the end of the list */ +}; + +#endif // KZGTEST \ No newline at end of file diff --git a/src/poly.h b/src/poly.h index 6d8b26a..ff0b5fe 100644 --- a/src/poly.h +++ b/src/poly.h @@ -33,11 +33,7 @@ typedef struct { } poly; void eval_poly(fr_t *out, const poly *p, const fr_t *x); -C_KZG_RET poly_long_div(poly *out, const poly *dividend, const poly *divisor); -C_KZG_RET poly_mul_direct(poly *out, const poly *a, const poly *b); -C_KZG_RET poly_mul_fft(poly *out, const poly *a, const poly *b, FFTSettings *fs); C_KZG_RET poly_inverse(poly *out, poly *b); -C_KZG_RET poly_fast_div(poly *out, const poly *dividend, const poly *divisor); C_KZG_RET poly_mul(poly *out, const poly *a, const poly *b); C_KZG_RET poly_mul_(poly *out, const poly *a, const poly *b, FFTSettings *fs); C_KZG_RET new_poly_div(poly *out, const poly *dividend, const poly *divisor); diff --git a/src/poly_test.c b/src/poly_test.c deleted file mode 100644 index 91dd1c1..0000000 --- a/src/poly_test.c +++ /dev/null @@ -1,416 +0,0 @@ -/* - * Copyright 2021 Benjamin Edgington - * - * Licensed under the Apache License, Version 2.0 (the "License"); - * you may not use this file except in compliance with the License. - * You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -#include "../inc/acutest.h" -#include "test_util.h" -#include "poly.h" - -// If anyone has a nicer way to initialise these test data, I'd love to hear it. -typedef struct polydata { - int length; - int coeffs[]; -} polydata; - -// (x^2 - 1) / (x + 1) = x - 1 -polydata test_0_0 = {3, {-1, 0, 1}}; -polydata test_0_1 = {2, {1, 1}}; -polydata test_0_2 = {2, {-1, 1}}; - -// (12x^3 - 11x^2 + 9x + 18) / (4x + 3) = 3x^2 - 5x + 6 -polydata test_1_0 = {4, {18, 9, -11, 12}}; -polydata test_1_1 = {2, {3, 4}}; -polydata test_1_2 = {3, {6, -5, 3}}; - -// (x + 1) / (x^2 - 1) = nil -polydata test_2_0 = {2, {1, 1}}; -polydata test_2_1 = {3, {-1, 0, 2}}; -polydata test_2_2 = {0, {}}; - -// (10x^2 + 20x + 30) / 10 = x^2 + 2x + 3 -polydata test_3_0 = {3, {30, 20, 10}}; -polydata test_3_1 = {1, {10}}; -polydata test_3_2 = {3, {3, 2, 1}}; - -// (x^2 + x) / (x + 1) = x -polydata test_4_0 = {3, {0, 1, 1}}; -polydata test_4_1 = {2, {1, 1}}; -polydata test_4_2 = {2, {0, 1}}; - -// (x^2 + x + 1) / 1 = x^2 + x + 1 -polydata test_5_0 = {3, {1, 1, 1}}; -polydata test_5_1 = {1, {1}}; -polydata test_5_2 = {3, {1, 1, 1}}; - -// (x^2 + x + 1) / (0x + 1) = x^2 + x + 1 -polydata test_6_0 = {3, {1, 1, 1}}; -polydata test_6_1 = {2, {1, 0}}; // The highest coefficient is zero -polydata test_6_2 = {3, {1, 1, 1}}; - -polydata *test[][3] = {{&test_0_0, &test_0_1, &test_0_2}, {&test_1_0, &test_1_1, &test_1_2}, - {&test_2_0, &test_2_1, &test_2_2}, {&test_3_0, &test_3_1, &test_3_2}, - {&test_4_0, &test_4_1, &test_4_2}, {&test_5_0, &test_5_1, &test_5_2}, - {&test_6_0, &test_6_1, &test_6_2}}; - -/* Internal utility function */ -void new_test_poly(poly *p, polydata *data) { - new_poly(p, data->length); - for (int i = 0; i < p->length; i++) { - int coeff = data->coeffs[i]; - if (coeff >= 0) { - fr_from_uint64(&p->coeffs[i], coeff); - } else { - fr_from_uint64(&p->coeffs[i], -coeff); - fr_negate(&p->coeffs[i], &p->coeffs[i]); - } - } -} - -void poly_test_div(void) { - poly dividend, divisor, expected, actual; - int ntest = sizeof test / sizeof test[0]; - - for (int i = 0; i < ntest; i++) { - new_test_poly(÷nd, test[i][0]); - new_test_poly(&divisor, test[i][1]); - new_test_poly(&expected, test[i][2]); - - if (TEST_CHECK(C_KZG_OK == new_poly_div(&actual, ÷nd, &divisor))) { - if (TEST_CHECK(actual.length == expected.length)) { - for (int j = 0; j < actual.length; j++) { - TEST_CHECK(fr_equal(&actual.coeffs[j], &expected.coeffs[j])); - TEST_MSG("Failed test %d with incorrect value", i); - } - } else { - TEST_MSG("Failed test %d with incorrect length.", i); - } - } else { - TEST_MSG("Failed test %d with bad return value.", i); - } - - free_poly(÷nd); - free_poly(&divisor); - free_poly(&expected); - free_poly(&actual); - } -} - -void poly_div_by_zero(void) { - fr_t a[2]; - poly dividend, divisor, dummy; - - // Calculate (x + 1) / 0 = FAIL - - // Dividend - fr_from_uint64(&a[0], 1); - fr_from_uint64(&a[1], 1); - dividend.length = 2; - dividend.coeffs = a; - - // Divisor - new_poly(&divisor, 0); - - TEST_CHECK(C_KZG_BADARGS == new_poly_div(&dummy, ÷nd, &divisor)); - - free_poly(&divisor); - free_poly(&dummy); -} - -void poly_eval_check(void) { - uint64_t n = 10; - fr_t actual, expected; - poly p; - new_poly(&p, n); - for (uint64_t i = 0; i < n; i++) { - fr_from_uint64(&p.coeffs[i], i + 1); - } - fr_from_uint64(&expected, n * (n + 1) / 2); - - eval_poly(&actual, &p, &fr_one); - - TEST_CHECK(fr_equal(&expected, &actual)); - - free_poly(&p); -} - -void poly_eval_0_check(void) { - uint64_t n = 7, a = 597; - fr_t actual, expected; - poly p; - new_poly(&p, n); - for (uint64_t i = 0; i < n; i++) { - fr_from_uint64(&p.coeffs[i], i + a); - } - fr_from_uint64(&expected, a); - - eval_poly(&actual, &p, &fr_zero); - - TEST_CHECK(fr_equal(&expected, &actual)); - - free_poly(&p); -} - -void poly_eval_nil_check(void) { - uint64_t n = 0; - fr_t actual; - poly p; - new_poly(&p, n); - - eval_poly(&actual, &p, &fr_one); - - TEST_CHECK(fr_equal(&fr_zero, &actual)); - - free_poly(&p); -} - -void poly_mul_direct_test(void) { - - // Calculate (3x^2 - 5x + 6) * (4x + 3) = 12x^3 - 11x^2 + 9x + 18 - static polydata multiplier_data = {3, {6, -5, 3}}; - static polydata multiplicand_data = {2, {3, 4}}; - static polydata expected_data = {4, {18, 9, -11, 12}}; - - poly multiplicand, multiplier, expected, actual0, actual1; - - new_test_poly(&multiplicand, &multiplicand_data); - new_test_poly(&multiplier, &multiplier_data); - new_test_poly(&expected, &expected_data); - - new_poly(&actual0, 4); - TEST_CHECK(C_KZG_OK == poly_mul_direct(&actual0, &multiplicand, &multiplier)); - TEST_CHECK(fr_equal(&expected.coeffs[0], &actual0.coeffs[0])); - TEST_CHECK(fr_equal(&expected.coeffs[1], &actual0.coeffs[1])); - TEST_CHECK(fr_equal(&expected.coeffs[2], &actual0.coeffs[2])); - TEST_CHECK(fr_equal(&expected.coeffs[3], &actual0.coeffs[3])); - - // Check commutativity - new_poly(&actual1, 4); - TEST_CHECK(C_KZG_OK == poly_mul_direct(&actual1, &multiplier, &multiplicand)); - TEST_CHECK(fr_equal(&expected.coeffs[0], &actual1.coeffs[0])); - TEST_CHECK(fr_equal(&expected.coeffs[1], &actual1.coeffs[1])); - TEST_CHECK(fr_equal(&expected.coeffs[2], &actual1.coeffs[2])); - TEST_CHECK(fr_equal(&expected.coeffs[3], &actual1.coeffs[3])); - - free_poly(&multiplicand); - free_poly(&multiplier); - free_poly(&expected); - free_poly(&actual0); - free_poly(&actual1); -} - -void poly_mul_fft_test(void) { - - // Calculate (3x^2 - 5x + 6) * (4x + 3) = 12x^3 - 11x^2 + 9x + 18 - static polydata multiplier_data = {3, {6, -5, 3}}; - static polydata multiplicand_data = {2, {3, 4}}; - static polydata expected_data = {4, {18, 9, -11, 12}}; - - poly multiplicand, multiplier, expected, actual0, actual1; - - new_test_poly(&multiplicand, &multiplicand_data); - new_test_poly(&multiplier, &multiplier_data); - new_test_poly(&expected, &expected_data); - - new_poly(&actual0, 4); - TEST_CHECK(C_KZG_OK == poly_mul_fft(&actual0, &multiplicand, &multiplier, NULL)); - TEST_CHECK(fr_equal(&expected.coeffs[0], &actual0.coeffs[0])); - TEST_CHECK(fr_equal(&expected.coeffs[1], &actual0.coeffs[1])); - TEST_CHECK(fr_equal(&expected.coeffs[2], &actual0.coeffs[2])); - TEST_CHECK(fr_equal(&expected.coeffs[3], &actual0.coeffs[3])); - - // Check commutativity - new_poly(&actual1, 4); - TEST_CHECK(C_KZG_OK == poly_mul_fft(&actual1, &multiplier, &multiplicand, NULL)); - TEST_CHECK(fr_equal(&expected.coeffs[0], &actual1.coeffs[0])); - TEST_CHECK(fr_equal(&expected.coeffs[1], &actual1.coeffs[1])); - TEST_CHECK(fr_equal(&expected.coeffs[2], &actual1.coeffs[2])); - TEST_CHECK(fr_equal(&expected.coeffs[3], &actual1.coeffs[3])); - - free_poly(&multiplicand); - free_poly(&multiplier); - free_poly(&expected); - free_poly(&actual0); - free_poly(&actual1); -} - -void poly_inverse_simple_0(void) { - - // 1 / (1 - x) = 1 + x + x^2 + ... - - poly p, q; - int d = 16; // The number of terms to take - - new_poly(&p, 2); - p.coeffs[0] = fr_one; - p.coeffs[1] = fr_one; - fr_negate(&p.coeffs[1], &p.coeffs[1]); - - new_poly(&q, d); - TEST_CHECK(C_KZG_OK == poly_inverse(&q, &p)); - - for (int i = 0; i < d; i++) { - TEST_CHECK(fr_is_one(&q.coeffs[i])); - } - - free_poly(&p); - free_poly(&q); -} - -void poly_inverse_simple_1(void) { - - // 1 / (1 + x) = 1 - x + x^2 - ... - - poly p, q; - int d = 16; // The number of terms to take - - new_poly(&p, 2); - p.coeffs[0] = fr_one; - p.coeffs[1] = fr_one; - - new_poly(&q, d); - TEST_CHECK(C_KZG_OK == poly_inverse(&q, &p)); - - for (int i = 0; i < d; i++) { - fr_t tmp = q.coeffs[i]; - if (i & 1) { - fr_negate(&tmp, &tmp); - } - TEST_CHECK(fr_is_one(&tmp)); - } - - free_poly(&p); - free_poly(&q); -} - -void poly_mul_random(void) { - - // Compare the output of poly_mul_direct() and poly_mul_fft() - - poly multiplicand, multiplier; - poly q0, q1; - - for (int k = 0; k < 256; k++) { - - int multiplicand_length = 1 + rand() % 1000; - int multiplier_length = 1 + rand() % 1000; - int out_length = 1 + rand() % 1000; - - new_poly(&multiplicand, multiplicand_length); - new_poly(&multiplier, multiplier_length); - - for (int i = 0; i < multiplicand_length; i++) { - multiplicand.coeffs[i] = rand_fr(); - } - for (int i = 0; i < multiplier_length; i++) { - multiplier.coeffs[i] = rand_fr(); - } - - // Ensure that the polynomials' orders corresponds to their lengths - if (fr_is_zero(&multiplicand.coeffs[multiplicand.length - 1])) { - multiplicand.coeffs[multiplicand.length - 1] = fr_one; - } - if (fr_is_zero(&multiplier.coeffs[multiplier.length - 1])) { - multiplier.coeffs[multiplier.length - 1] = fr_one; - } - - new_poly(&q0, out_length); // Truncate the result - TEST_CHECK(C_KZG_OK == poly_mul_direct(&q0, &multiplicand, &multiplier)); - - new_poly(&q1, out_length); - TEST_CHECK(C_KZG_OK == poly_mul_fft(&q1, &multiplicand, &multiplier, NULL)); - - TEST_CHECK(q1.length == q0.length); - for (int i = 0; i < q0.length; i++) { - if (!TEST_CHECK(fr_equal(&q0.coeffs[i], &q1.coeffs[i]))) { - TEST_MSG("round = %d, i = %d, multiplicand_length = %lu, multiplier_length = %lu, out_length = %lu", k, - i, multiplicand.length, multiplier.length, q0.length); - } - } - - free_poly(&multiplicand); - free_poly(&multiplier); - free_poly(&q0); - free_poly(&q1); - } -} - -void poly_div_random(void) { - - // Compare the output of poly_fast_div() and poly_long_div() - - poly dividend, divisor; - poly q0, q1; - - for (int k = 0; k < 256; k++) { - - int dividend_length = 2 + rand() % 1000; - int divisor_length = 1 + rand() % dividend_length; - - new_poly(÷nd, dividend_length); - new_poly(&divisor, divisor_length); - - for (int i = 0; i < dividend_length; i++) { - dividend.coeffs[i] = rand_fr(); - } - for (int i = 0; i < divisor_length; i++) { - divisor.coeffs[i] = rand_fr(); - } - - // Ensure that the polynomials' orders corresponds to their lengths - if (fr_is_zero(÷nd.coeffs[dividend.length - 1])) { - dividend.coeffs[dividend.length - 1] = fr_one; - } - if (fr_is_zero(&divisor.coeffs[divisor.length - 1])) { - divisor.coeffs[divisor.length - 1] = fr_one; - } - - new_poly(&q0, dividend.length - divisor.length + 1); - TEST_CHECK(C_KZG_OK == poly_long_div(&q0, ÷nd, &divisor)); - - new_poly(&q1, dividend.length - divisor.length + 1); - TEST_CHECK(C_KZG_OK == poly_fast_div(&q1, ÷nd, &divisor)); - - TEST_CHECK(q1.length == q0.length); - for (int i = 0; i < q0.length; i++) { - if (!TEST_CHECK(fr_equal(&q0.coeffs[i], &q1.coeffs[i]))) { - TEST_MSG("round = %d, dividend_length = %lu, divisor_length = %lu, i = %d", k, dividend.length, - divisor.length, i); - } - } - - free_poly(÷nd); - free_poly(&divisor); - free_poly(&q0); - free_poly(&q1); - } -} - -TEST_LIST = { - {"POLY_TEST", title}, - {"poly_test_div", poly_test_div}, -#ifndef DEBUG - {"poly_div_by_zero", poly_div_by_zero}, -#endif - {"poly_eval_check", poly_eval_check}, - {"poly_eval_0_check", poly_eval_0_check}, - {"poly_eval_nil_check", poly_eval_nil_check}, - {"poly_mul_direct_test", poly_mul_direct_test}, - {"poly_mul_fft_test", poly_mul_fft_test}, - {"poly_inverse_simple_0", poly_inverse_simple_0}, - {"poly_inverse_simple_1", poly_inverse_simple_1}, - {"poly_mul_random", poly_mul_random}, - {"poly_div_random", poly_div_random}, - {NULL, NULL} /* zero record marks the end of the list */ -}; diff --git a/src/recover.c b/src/recover.c index 920a4a5..691ada3 100644 --- a/src/recover.c +++ b/src/recover.c @@ -38,7 +38,7 @@ * @param[out,in] p The polynomial coefficients to be scaled * @param[in] len_p Length of the polynomial coefficients */ -void scale_poly(fr_t *p, uint64_t len_p) { +static void scale_poly(fr_t *p, uint64_t len_p) { fr_t scale_factor, factor_power, inv_factor; fr_from_uint64(&scale_factor, SCALE_FACTOR); fr_inv(&inv_factor, &scale_factor); @@ -59,7 +59,7 @@ void scale_poly(fr_t *p, uint64_t len_p) { * @param[out,in] p The polynomial coefficients to be unscaled * @param[in] len_p Length of the polynomial coefficients */ -void unscale_poly(fr_t *p, uint64_t len_p) { +static void unscale_poly(fr_t *p, uint64_t len_p) { fr_t scale_factor, factor_power; fr_from_uint64(&scale_factor, SCALE_FACTOR); factor_power = fr_one; @@ -177,4 +177,130 @@ C_KZG_RET recover_poly_from_samples(fr_t *reconstructed_data, fr_t *samples, uin free(missing); return C_KZG_OK; -} \ No newline at end of file +} + +#ifdef KZGTEST + +#include "../inc/acutest.h" +#include "test_util.h" + +// Utility for setting a random (len_data - known) number of elements to NULL +void random_missing(fr_t *with_missing, fr_t *data, uint64_t len_data, uint64_t known) { + uint64_t *missing_idx; + TEST_CHECK(C_KZG_OK == new_uint64_array(&missing_idx, len_data)); + for (uint64_t i = 0; i < len_data; i++) { + missing_idx[i] = i; + } + shuffle(missing_idx, len_data); + + for (uint64_t i = 0; i < len_data; i++) { + with_missing[i] = data[i]; + } + for (uint64_t i = 0; i < len_data - known; i++) { + with_missing[missing_idx[i]] = fr_null; + } + + free(missing_idx); +} + +void recover_simple(void) { + FFTSettings fs; + TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 2)); + + fr_t poly[fs.max_width]; + for (int i = 0; i < fs.max_width / 2; i++) { + fr_from_uint64(&poly[i], i); + } + for (int i = fs.max_width / 2; i < fs.max_width; i++) { + poly[i] = fr_zero; + } + + fr_t data[fs.max_width]; + TEST_CHECK(C_KZG_OK == fft_fr(data, poly, false, fs.max_width, &fs)); + + fr_t sample[fs.max_width]; + sample[0] = data[0]; + sample[1] = fr_null; + sample[2] = fr_null; + sample[3] = data[3]; + + fr_t recovered[fs.max_width]; + TEST_CHECK(C_KZG_OK == recover_poly_from_samples(recovered, sample, fs.max_width, &fs)); + + // Check recovered data + for (int i = 0; i < fs.max_width; i++) { + TEST_CHECK(fr_equal(&data[i], &recovered[i])); + } + + // Also check against original coefficients + fr_t back[fs.max_width]; + TEST_CHECK(C_KZG_OK == fft_fr(back, recovered, true, fs.max_width, &fs)); + for (int i = 0; i < fs.max_width / 2; i++) { + TEST_CHECK(fr_equal(&poly[i], &back[i])); + } + for (int i = fs.max_width / 2; i < fs.max_width; i++) { + TEST_CHECK(fr_is_zero(&back[i])); + } + + free_fft_settings(&fs); +} + +void recover_random(void) { + FFTSettings fs; + TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 12)); + + fr_t *poly, *data, *samples, *recovered, *back; + TEST_CHECK(C_KZG_OK == new_fr_array(&poly, fs.max_width)); + TEST_CHECK(C_KZG_OK == new_fr_array(&data, fs.max_width)); + TEST_CHECK(C_KZG_OK == new_fr_array(&samples, fs.max_width)); + TEST_CHECK(C_KZG_OK == new_fr_array(&recovered, fs.max_width)); + TEST_CHECK(C_KZG_OK == new_fr_array(&back, fs.max_width)); + + for (int i = 0; i < fs.max_width / 2; i++) { + fr_from_uint64(&poly[i], i); + } + for (int i = fs.max_width / 2; i < fs.max_width; i++) { + poly[i] = fr_zero; + } + + TEST_CHECK(C_KZG_OK == fft_fr(data, poly, false, fs.max_width, &fs)); + + // Having half of the data is the minimum + for (float known_ratio = 0.5; known_ratio < 1.0; known_ratio += 0.05) { + uint64_t known = fs.max_width * known_ratio; + for (int i = 0; i < 4; i++) { + random_missing(samples, data, fs.max_width, known); + + TEST_CHECK(C_KZG_OK == recover_poly_from_samples(recovered, samples, fs.max_width, &fs)); + for (int i = 0; i < fs.max_width; i++) { + TEST_CHECK(fr_equal(&data[i], &recovered[i])); + } + + // Also check against original coefficients + fr_t back[fs.max_width]; + TEST_CHECK(C_KZG_OK == fft_fr(back, recovered, true, fs.max_width, &fs)); + for (int i = 0; i < fs.max_width / 2; i++) { + TEST_CHECK(fr_equal(&poly[i], &back[i])); + } + for (int i = fs.max_width / 2; i < fs.max_width; i++) { + TEST_CHECK(fr_is_zero(&back[i])); + } + } + } + + free(poly); + free(data); + free(samples); + free(recovered); + free(back); + free_fft_settings(&fs); +} + +TEST_LIST = { + {"RECOVER_TEST", title}, + {"recover_simple", recover_simple}, + {"recover_random", recover_random}, + {NULL, NULL} /* zero record marks the end of the list */ +}; + +#endif // KZGTEST \ No newline at end of file diff --git a/src/recover_test.c b/src/recover_test.c deleted file mode 100644 index 91a61b3..0000000 --- a/src/recover_test.c +++ /dev/null @@ -1,141 +0,0 @@ -/* - * Copyright 2021 Benjamin Edgington - * - * Licensed under the Apache License, Version 2.0 (the "License"); - * you may not use this file except in compliance with the License. - * You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -#include "../inc/acutest.h" -#include "c_kzg_util.h" -#include "test_util.h" -#include "recover.h" -#include "fft_fr.h" -#include "debug_util.h" - -// Utility for setting a random (len_data - known) number of elements to NULL -void random_missing(fr_t *with_missing, fr_t *data, uint64_t len_data, uint64_t known) { - uint64_t *missing_idx; - TEST_CHECK(C_KZG_OK == new_uint64_array(&missing_idx, len_data)); - for (uint64_t i = 0; i < len_data; i++) { - missing_idx[i] = i; - } - shuffle(missing_idx, len_data); - - for (uint64_t i = 0; i < len_data; i++) { - with_missing[i] = data[i]; - } - for (uint64_t i = 0; i < len_data - known; i++) { - with_missing[missing_idx[i]] = fr_null; - } - - free(missing_idx); -} - -void recover_simple(void) { - FFTSettings fs; - TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 2)); - - fr_t poly[fs.max_width]; - for (int i = 0; i < fs.max_width / 2; i++) { - fr_from_uint64(&poly[i], i); - } - for (int i = fs.max_width / 2; i < fs.max_width; i++) { - poly[i] = fr_zero; - } - - fr_t data[fs.max_width]; - TEST_CHECK(C_KZG_OK == fft_fr(data, poly, false, fs.max_width, &fs)); - - fr_t sample[fs.max_width]; - sample[0] = data[0]; - sample[1] = fr_null; - sample[2] = fr_null; - sample[3] = data[3]; - - fr_t recovered[fs.max_width]; - TEST_CHECK(C_KZG_OK == recover_poly_from_samples(recovered, sample, fs.max_width, &fs)); - - // Check recovered data - for (int i = 0; i < fs.max_width; i++) { - TEST_CHECK(fr_equal(&data[i], &recovered[i])); - } - - // Also check against original coefficients - fr_t back[fs.max_width]; - TEST_CHECK(C_KZG_OK == fft_fr(back, recovered, true, fs.max_width, &fs)); - for (int i = 0; i < fs.max_width / 2; i++) { - TEST_CHECK(fr_equal(&poly[i], &back[i])); - } - for (int i = fs.max_width / 2; i < fs.max_width; i++) { - TEST_CHECK(fr_is_zero(&back[i])); - } - - free_fft_settings(&fs); -} - -void recover_random(void) { - FFTSettings fs; - TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 12)); - - fr_t *poly, *data, *samples, *recovered, *back; - TEST_CHECK(C_KZG_OK == new_fr_array(&poly, fs.max_width)); - TEST_CHECK(C_KZG_OK == new_fr_array(&data, fs.max_width)); - TEST_CHECK(C_KZG_OK == new_fr_array(&samples, fs.max_width)); - TEST_CHECK(C_KZG_OK == new_fr_array(&recovered, fs.max_width)); - TEST_CHECK(C_KZG_OK == new_fr_array(&back, fs.max_width)); - - for (int i = 0; i < fs.max_width / 2; i++) { - fr_from_uint64(&poly[i], i); - } - for (int i = fs.max_width / 2; i < fs.max_width; i++) { - poly[i] = fr_zero; - } - - TEST_CHECK(C_KZG_OK == fft_fr(data, poly, false, fs.max_width, &fs)); - - // Having half of the data is the minimum - for (float known_ratio = 0.5; known_ratio < 1.0; known_ratio += 0.05) { - uint64_t known = fs.max_width * known_ratio; - for (int i = 0; i < 4; i++) { - random_missing(samples, data, fs.max_width, known); - - TEST_CHECK(C_KZG_OK == recover_poly_from_samples(recovered, samples, fs.max_width, &fs)); - for (int i = 0; i < fs.max_width; i++) { - TEST_CHECK(fr_equal(&data[i], &recovered[i])); - } - - // Also check against original coefficients - fr_t back[fs.max_width]; - TEST_CHECK(C_KZG_OK == fft_fr(back, recovered, true, fs.max_width, &fs)); - for (int i = 0; i < fs.max_width / 2; i++) { - TEST_CHECK(fr_equal(&poly[i], &back[i])); - } - for (int i = fs.max_width / 2; i < fs.max_width; i++) { - TEST_CHECK(fr_is_zero(&back[i])); - } - } - } - - free(poly); - free(data); - free(samples); - free(recovered); - free(back); - free_fft_settings(&fs); -} - -TEST_LIST = { - {"RECOVER_TEST", title}, - {"recover_simple", recover_simple}, - {"recover_random", recover_random}, - {NULL, NULL} /* zero record marks the end of the list */ -}; diff --git a/src/utility.c b/src/utility.c index 64a5c2a..212f525 100644 --- a/src/utility.c +++ b/src/utility.c @@ -152,3 +152,173 @@ C_KZG_RET reverse_bit_order(void *values, size_t size, uint64_t n) { return C_KZG_OK; } + +#ifdef KZGTEST + +#include "../inc/acutest.h" +#include "test_util.h" +#include "c_kzg_util.h" + +static uint32_t rev_bits_slow(uint32_t a) { + uint32_t ret = 0; + for (int i = 0; i < 32; i++) { + ret <<= 1; + ret |= a & 1; + a >>= 1; + } + return ret; +} + +void is_power_of_two_works(void) { + // All actual powers of two + for (int i = 0; i <= 63; i++) { + TEST_CHECK(true == is_power_of_two((uint64_t)1 << i)); + TEST_MSG("Case %d", i); + } + + // This is a bit weird + TEST_CHECK(true == is_power_of_two(0)); + + // Not powers of two + TEST_CHECK(false == is_power_of_two(123)); + TEST_CHECK(false == is_power_of_two(1234567)); +} + +void test_log2_pow2(void) { + int actual, expected; + for (int i = 0; i < 32; i++) { + expected = i; + actual = log2_pow2((uint32_t)1 << i); + TEST_CHECK(expected == actual); + } +} + +void test_next_power_of_two_powers(void) { + for (int i = 0; i <= 63; i++) { + uint64_t expected = (uint64_t)1 << i; + uint64_t actual = next_power_of_two(expected); + TEST_CHECK(expected == actual); + } +} + +void test_next_power_of_two_random(void) { + for (int i = 0; i < 32768; i++) { + uint64_t a = 1 + (rand_uint64() >> 1); // It's not expected to work for a > 2^63 + uint64_t higher = next_power_of_two(a); + uint64_t lower = higher >> 1; + if (!(TEST_CHECK(is_power_of_two(higher)) && TEST_CHECK(higher >= a) && TEST_CHECK(lower < a))) { + TEST_MSG("Failed for %lu", a); + } + } +} + +void test_reverse_bits_macros(void) { + TEST_CHECK(128 == rev_byte(1)); + TEST_CHECK(128 == rev_byte(257)); + TEST_CHECK((uint32_t)1 << 31 == rev_4byte(1)); + TEST_CHECK(0x1e6a2c48 == rev_4byte(0x12345678)); + TEST_CHECK(0x00000000 == rev_4byte(0x00000000)); + TEST_CHECK(0xffffffff == rev_4byte(0xffffffff)); +} + +void test_reverse_bits_powers(void) { + uint32_t actual, expected; + for (int i = 0; i < 32; i++) { + expected = (uint32_t)1 << (31 - i); + actual = reverse_bits((uint32_t)1 << i); + TEST_CHECK(expected == actual); + } +} + +void test_reverse_bits_random(void) { + for (int i = 0; i < 32768; i++) { + uint32_t a = (uint32_t)(rand_uint64() & 0xffffffffL); + TEST_CHECK(rev_bits_slow(a) == reverse_bits(a)); + TEST_MSG("Failed for %08x. Expected %08x, got %08x.", a, rev_bits_slow(a), reverse_bits(a)); + } +} + +void test_reverse_bit_order_g1(void) { + int size = 10, n = 1 << size; + g1_t a[n], b[n]; + fr_t tmp; + + for (int i = 0; i < n; i++) { + fr_from_uint64(&tmp, i); + g1_mul(&a[i], &g1_generator, &tmp); + b[i] = a[i]; + } + + TEST_CHECK(C_KZG_OK == reverse_bit_order(a, sizeof(g1_t), n)); + for (int i = 0; i < n; i++) { + TEST_CHECK(true == g1_equal(&b[reverse_bits(i) >> (32 - size)], &a[i])); + } + + // Hand check a few select values + TEST_CHECK(true == g1_equal(&b[0], &a[0])); + TEST_CHECK(false == g1_equal(&b[1], &a[1])); + TEST_CHECK(true == g1_equal(&b[n - 1], &a[n - 1])); +} + +void test_reverse_bit_order_fr(void) { + int size = 12, n = 1 << size; + fr_t a[n], b[n]; + + for (int i = 0; i < n; i++) { + fr_from_uint64(&a[i], i); + b[i] = a[i]; + } + + TEST_CHECK(C_KZG_OK == reverse_bit_order(a, sizeof(fr_t), n)); + for (int i = 0; i < n; i++) { + TEST_CHECK(true == fr_equal(&b[reverse_bits(i) >> (32 - size)], &a[i])); + } + + // Hand check a few select values + TEST_CHECK(true == fr_equal(&b[0], &a[0])); + TEST_CHECK(false == fr_equal(&b[1], &a[1])); + TEST_CHECK(true == fr_equal(&b[n - 1], &a[n - 1])); +} + +void test_reverse_bit_order_fr_large(void) { + int size = 22, n = 1 << size; + fr_t *a, *b; + + TEST_CHECK(C_KZG_OK == new_fr_array(&a, n)); + TEST_CHECK(C_KZG_OK == new_fr_array(&b, n)); + + for (int i = 0; i < n; i++) { + fr_from_uint64(&a[i], i); + b[i] = a[i]; + } + + TEST_CHECK(C_KZG_OK == reverse_bit_order(a, sizeof(fr_t), n)); + for (int i = 0; i < n; i++) { + TEST_CHECK(true == fr_equal(&b[reverse_bits(i) >> (32 - size)], &a[i])); + } + + // Hand check a few select values + TEST_CHECK(true == fr_equal(&b[0], &a[0])); + TEST_CHECK(false == fr_equal(&b[1], &a[1])); + TEST_CHECK(true == fr_equal(&b[n - 1], &a[n - 1])); + + free(a); + free(b); +} + +TEST_LIST = { + {"UTILITY_TEST", title}, + {"is_power_of_two_works", is_power_of_two_works}, + {"test_log2_pow2", test_log2_pow2}, + {"test_next_power_of_two_powers", test_next_power_of_two_powers}, + {"test_next_power_of_two_random", test_next_power_of_two_random}, + {"test_reverse_bits_macros", test_reverse_bits_macros}, + {"test_reverse_bits_powers", test_reverse_bits_powers}, + {"test_reverse_bits_random", test_reverse_bits_random}, + {"test_reverse_bit_order_g1", test_reverse_bit_order_g1}, + {"test_reverse_bit_order_fr", test_reverse_bit_order_fr}, + {"test_reverse_bit_order_fr_large", test_reverse_bit_order_fr_large}, + {NULL, NULL} /* zero record marks the end of the list */ +}; + +#endif // KZGTEST \ No newline at end of file diff --git a/src/utility_test.c b/src/utility_test.c deleted file mode 100644 index 9f40df9..0000000 --- a/src/utility_test.c +++ /dev/null @@ -1,182 +0,0 @@ -/* - * Copyright 2021 Benjamin Edgington - * - * Licensed under the Apache License, Version 2.0 (the "License"); - * you may not use this file except in compliance with the License. - * You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -#include "../inc/acutest.h" -#include "test_util.h" -#include "utility.h" -#include "c_kzg_util.h" - -static uint32_t rev_bits_slow(uint32_t a) { - uint32_t ret = 0; - for (int i = 0; i < 32; i++) { - ret <<= 1; - ret |= a & 1; - a >>= 1; - } - return ret; -} - -void is_power_of_two_works(void) { - // All actual powers of two - for (int i = 0; i <= 63; i++) { - TEST_CHECK(true == is_power_of_two((uint64_t)1 << i)); - TEST_MSG("Case %d", i); - } - - // This is a bit weird - TEST_CHECK(true == is_power_of_two(0)); - - // Not powers of two - TEST_CHECK(false == is_power_of_two(123)); - TEST_CHECK(false == is_power_of_two(1234567)); -} - -void test_log2_pow2(void) { - int actual, expected; - for (int i = 0; i < 32; i++) { - expected = i; - actual = log2_pow2((uint32_t)1 << i); - TEST_CHECK(expected == actual); - } -} - -void test_next_power_of_two_powers(void) { - for (int i = 0; i <= 63; i++) { - uint64_t expected = (uint64_t)1 << i; - uint64_t actual = next_power_of_two(expected); - TEST_CHECK(expected == actual); - } -} - -void test_next_power_of_two_random(void) { - for (int i = 0; i < 32768; i++) { - uint64_t a = 1 + (rand_uint64() >> 1); // It's not expected to work for a > 2^63 - uint64_t higher = next_power_of_two(a); - uint64_t lower = higher >> 1; - if (!(TEST_CHECK(is_power_of_two(higher)) && TEST_CHECK(higher >= a) && TEST_CHECK(lower < a))) { - TEST_MSG("Failed for %lu", a); - } - } -} - -void test_reverse_bits_macros(void) { - TEST_CHECK(128 == rev_byte(1)); - TEST_CHECK(128 == rev_byte(257)); - TEST_CHECK((uint32_t)1 << 31 == rev_4byte(1)); - TEST_CHECK(0x1e6a2c48 == rev_4byte(0x12345678)); - TEST_CHECK(0x00000000 == rev_4byte(0x00000000)); - TEST_CHECK(0xffffffff == rev_4byte(0xffffffff)); -} - -void test_reverse_bits_powers(void) { - uint32_t actual, expected; - for (int i = 0; i < 32; i++) { - expected = (uint32_t)1 << (31 - i); - actual = reverse_bits((uint32_t)1 << i); - TEST_CHECK(expected == actual); - } -} - -void test_reverse_bits_random(void) { - for (int i = 0; i < 32768; i++) { - uint32_t a = (uint32_t)(rand_uint64() & 0xffffffffL); - TEST_CHECK(rev_bits_slow(a) == reverse_bits(a)); - TEST_MSG("Failed for %08x. Expected %08x, got %08x.", a, rev_bits_slow(a), reverse_bits(a)); - } -} - -void test_reverse_bit_order_g1(void) { - int size = 10, n = 1 << size; - g1_t a[n], b[n]; - fr_t tmp; - - for (int i = 0; i < n; i++) { - fr_from_uint64(&tmp, i); - g1_mul(&a[i], &g1_generator, &tmp); - b[i] = a[i]; - } - - TEST_CHECK(C_KZG_OK == reverse_bit_order(a, sizeof(g1_t), n)); - for (int i = 0; i < n; i++) { - TEST_CHECK(true == g1_equal(&b[reverse_bits(i) >> (32 - size)], &a[i])); - } - - // Hand check a few select values - TEST_CHECK(true == g1_equal(&b[0], &a[0])); - TEST_CHECK(false == g1_equal(&b[1], &a[1])); - TEST_CHECK(true == g1_equal(&b[n - 1], &a[n - 1])); -} - -void test_reverse_bit_order_fr(void) { - int size = 12, n = 1 << size; - fr_t a[n], b[n]; - - for (int i = 0; i < n; i++) { - fr_from_uint64(&a[i], i); - b[i] = a[i]; - } - - TEST_CHECK(C_KZG_OK == reverse_bit_order(a, sizeof(fr_t), n)); - for (int i = 0; i < n; i++) { - TEST_CHECK(true == fr_equal(&b[reverse_bits(i) >> (32 - size)], &a[i])); - } - - // Hand check a few select values - TEST_CHECK(true == fr_equal(&b[0], &a[0])); - TEST_CHECK(false == fr_equal(&b[1], &a[1])); - TEST_CHECK(true == fr_equal(&b[n - 1], &a[n - 1])); -} - -void test_reverse_bit_order_fr_large(void) { - int size = 22, n = 1 << size; - fr_t *a, *b; - - TEST_CHECK(C_KZG_OK == new_fr_array(&a, n)); - TEST_CHECK(C_KZG_OK == new_fr_array(&b, n)); - - for (int i = 0; i < n; i++) { - fr_from_uint64(&a[i], i); - b[i] = a[i]; - } - - TEST_CHECK(C_KZG_OK == reverse_bit_order(a, sizeof(fr_t), n)); - for (int i = 0; i < n; i++) { - TEST_CHECK(true == fr_equal(&b[reverse_bits(i) >> (32 - size)], &a[i])); - } - - // Hand check a few select values - TEST_CHECK(true == fr_equal(&b[0], &a[0])); - TEST_CHECK(false == fr_equal(&b[1], &a[1])); - TEST_CHECK(true == fr_equal(&b[n - 1], &a[n - 1])); - - free(a); - free(b); -} - -TEST_LIST = { - {"UTILITY_TEST", title}, - {"is_power_of_two_works", is_power_of_two_works}, - {"test_log2_pow2", test_log2_pow2}, - {"test_next_power_of_two_powers", test_next_power_of_two_powers}, - {"test_next_power_of_two_random", test_next_power_of_two_random}, - {"test_reverse_bits_macros", test_reverse_bits_macros}, - {"test_reverse_bits_powers", test_reverse_bits_powers}, - {"test_reverse_bits_random", test_reverse_bits_random}, - {"test_reverse_bit_order_g1", test_reverse_bit_order_g1}, - {"test_reverse_bit_order_fr", test_reverse_bit_order_fr}, - {"test_reverse_bit_order_fr_large", test_reverse_bit_order_fr_large}, - {NULL, NULL} /* zero record marks the end of the list */ -}; diff --git a/src/zero_poly.c b/src/zero_poly.c index cdd7cd0..ade6f28 100644 --- a/src/zero_poly.c +++ b/src/zero_poly.c @@ -40,8 +40,8 @@ * @retval C_CZK_OK All is well * @retval C_CZK_BADARGS Invalid parameters were supplied */ -C_KZG_RET do_zero_poly_mul_partial(poly *dst, const uint64_t *indices, uint64_t len_indices, uint64_t stride, - const FFTSettings *fs) { +static C_KZG_RET do_zero_poly_mul_partial(poly *dst, const uint64_t *indices, uint64_t len_indices, uint64_t stride, + const FFTSettings *fs) { CHECK(len_indices > 0); CHECK(dst->length >= len_indices + 1); @@ -80,7 +80,7 @@ C_KZG_RET do_zero_poly_mul_partial(poly *dst, const uint64_t *indices, uint64_t * @retval C_CZK_OK All is well * @retval C_CZK_BADARGS Invalid parameters were supplied */ -C_KZG_RET pad_p(fr_t *out, uint64_t out_len, const poly *p) { +static C_KZG_RET pad_p(fr_t *out, uint64_t out_len, const poly *p) { CHECK(out_len >= p->length); for (uint64_t i = 0; i < p->length; i++) { out[i] = p->coeffs[i]; @@ -108,8 +108,8 @@ C_KZG_RET pad_p(fr_t *out, uint64_t out_len, const poly *p) { * @retval C_CZK_BADARGS Invalid parameters were supplied * @retval C_CZK_ERROR An internal error occurred */ -C_KZG_RET reduce_partials(poly *out, uint64_t len_out, fr_t *scratch, uint64_t len_scratch, const poly *partials, - uint64_t partial_count, const FFTSettings *fs) { +static C_KZG_RET reduce_partials(poly *out, uint64_t len_out, fr_t *scratch, uint64_t len_scratch, const poly *partials, + uint64_t partial_count, const FFTSettings *fs) { CHECK(is_power_of_two(len_out)); CHECK(len_scratch >= 3 * len_out); CHECK(partial_count > 0); @@ -255,3 +255,402 @@ C_KZG_RET zero_polynomial_via_multiplication(fr_t *zero_eval, poly *zero_poly, u return C_KZG_OK; } + +#ifdef KZGTEST + +#include "../inc/acutest.h" +#include "test_util.h" + +bool exists[16] = {true, false, false, true, false, true, true, false, + false, false, true, true, false, true, false, true}; + +uint64_t expected_eval_u64[16][4] = { + {0xfd5a5130b97ce0c3L, 0xb4748a4cb0f90e6dL, 0x12a1ab34b25b18c1L, 0x5a5ac0c81c9f7ea8L}, + {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, + {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, + {0xaa385cbce3dd1657L, 0x2fdab57a38bdb514L, 0x20e022e205dafa53L, 0x14077dd3f5d996b1L}, + {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, + {0x194018614b6f7276L, 0xdf2b18f870532376L, 0x1ff427cd5b583fe6L, 0x014d6444ff03dd09L}, + {0xcc84c2de684c0ddeL, 0xf1e7ab32aa830d02L, 0x967bf35a2a691f20L, 0x046109731cdf0d3cL}, + {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, + {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, + {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, + {0x96cddd2924212afbL, 0xeaa4c1f51421d8d8L, 0x3ae969cfa34d0ed1L, 0x6b6c5e876bc3916dL}, + {0x449310802f74ad49L, 0x47c940979163037aL, 0x10d311564afb9b2aL, 0x269b8531c369bafbL}, + {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, + {0xd9af75fe35c16cf1L, 0x068bb140cea92f75L, 0xe769811965e10a47L, 0x48ed97e6745612f2L}, + {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, + {0x7ef1f59bb1677685L, 0x33a637296680e8ceL, 0xaaf62b3f6e016709L, 0x454a299178a4dba9L}}; + +uint64_t expected_poly_u64[16][4] = { + {0xac159e2688bd4333L, 0x3bfef0f00df2ec88L, 0x561dcd0fd4d314d9L, 0x533bd8c1e977024eL}, + {0x18bc6eedc010ef8dL, 0xc731a3eb4ea2ab70L, 0x5c2589357ae121a8L, 0x04f9108d308f7016L}, + {0x232759f49556ac08L, 0x9776fe2e9f4c613cL, 0x74d5bed4eb2de960L, 0x1f6cf6719bfa0e68L}, + {0xf2f3461e8ab1ae34L, 0xeb220fcc11ef1c80L, 0x7a4637d3a637739bL, 0x19901a58cd177c53L}, + {0x9340f62465a1f4feL, 0xd9cb3ea6de494a11L, 0xee92ebc763cdff5dL, 0x5443e89811b5b9f5L}, + {0x269a255e2e4e48a4L, 0xfadae7a89d9b2f2bL, 0xb5515799b41e1a88L, 0x2e990979a0ffcee5L}, + {0x1c2f3a5759088c29L, 0x2a958d654cf1795fL, 0x9ca121fa43d152d1L, 0x1425239535953093L}, + {0x4c634e2d63ad89fdL, 0xd6ea7bc7da4ebe1aL, 0x9730a8fb88c7c895L, 0x1a01ffae0477c2a8L}, + {0x0000000000000001L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, + {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, + {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, + {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, + {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, + {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, + {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, + {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}}; + +void test_reduce_partials(void) { + FFTSettings fs; + TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 4)); + fr_t from_tree_reduction_coeffs[16], from_direct_coeffs[9], scratch[48]; + poly from_tree_reduction, from_direct; + from_tree_reduction.coeffs = from_tree_reduction_coeffs; + from_tree_reduction.length = 16; + from_direct.coeffs = from_direct_coeffs; + from_direct.length = 9; + + // Via reduce_partials + + poly partials[4]; + fr_t partial0[3], partial1[3], partial2[3], partial3[3]; + partials[0].coeffs = partial0, partials[0].length = 3; + partials[1].coeffs = partial1, partials[1].length = 3; + partials[2].coeffs = partial2, partials[2].length = 3; + partials[3].coeffs = partial3, partials[3].length = 3; + const uint64_t partial_indices[4][2] = {{1, 3}, {7, 8}, {9, 10}, {12, 13}}; + for (int i = 0; i < 4; i++) { + TEST_CHECK(C_KZG_OK == do_zero_poly_mul_partial(&partials[i], partial_indices[i], 2, 1, &fs)); + } + TEST_CHECK(C_KZG_OK == reduce_partials(&from_tree_reduction, 16, scratch, 48, partials, 4, &fs)); + + // Direct + uint64_t indices[] = {1, 3, 7, 8, 9, 10, 12, 13}; + TEST_CHECK(C_KZG_OK == do_zero_poly_mul_partial(&from_direct, indices, 8, 1, &fs)); + + // Compare + for (int i = 0; i < 9; i++) { + TEST_CHECK(fr_equal(&from_tree_reduction.coeffs[i], &from_direct.coeffs[i])); + } + + free_fft_settings(&fs); +} + +void reduce_partials_random(void) { + for (int scale = 5; scale < 13; scale++) { + for (int ii = 1; ii <= 7; ii++) { + float missing_ratio = 0.1 * ii; + + FFTSettings fs; + new_fft_settings(&fs, scale); + uint64_t point_count = fs.max_width; + uint64_t missing_count = point_count * missing_ratio; + + uint64_t *missing; + TEST_CHECK(C_KZG_OK == new_uint64_array(&missing, point_count)); + for (uint64_t i = 0; i < point_count; i++) { + missing[i] = i; + } + shuffle(missing, point_count); + + // Build the partials + poly *partials; + const int missing_per_partial = 63; + uint64_t indices[missing_per_partial]; + uint64_t partial_count = (missing_count + missing_per_partial - 1) / missing_per_partial; + TEST_CHECK(C_KZG_OK == new_poly_array(&partials, partial_count)); + for (uint64_t i = 0; i < partial_count; i++) { + uint64_t start = i * missing_per_partial; + uint64_t end = start + missing_per_partial; + if (end > missing_count) end = missing_count; + uint64_t partial_size = end - start; + TEST_CHECK(C_KZG_OK == new_fr_array(&partials[i].coeffs, partial_size + 1)); + for (int j = 0; j < partial_size; j++) { + indices[j] = missing[i * missing_per_partial + j]; + } + partials[i].length = partial_size + 1; + TEST_CHECK(C_KZG_OK == do_zero_poly_mul_partial(&partials[i], indices, partial_size, 1, &fs)); + } + + // From tree reduction + poly from_tree_reduction; + TEST_CHECK(C_KZG_OK == new_poly(&from_tree_reduction, point_count)); + fr_t *scratch; + TEST_CHECK(C_KZG_OK == new_fr_array(&scratch, point_count * 3)); + TEST_CHECK(C_KZG_OK == reduce_partials(&from_tree_reduction, point_count, scratch, point_count * 3, + partials, partial_count, &fs)); + + // From direct + poly from_direct; + TEST_CHECK(C_KZG_OK == new_poly(&from_direct, missing_count + 1)); + TEST_CHECK(C_KZG_OK == + do_zero_poly_mul_partial(&from_direct, missing, missing_count, fs.max_width / point_count, &fs)); + + for (uint64_t i = 0; i < missing_count + 1; i++) { + TEST_CHECK(fr_equal(&from_tree_reduction.coeffs[i], &from_direct.coeffs[i])); + } + + free_poly(&from_tree_reduction); + free_poly(&from_direct); + free(scratch); + for (uint64_t i = 0; i < partial_count; i++) { + free_poly(&partials[i]); + } + free(partials); + free(missing); + free_fft_settings(&fs); + } + } +} + +void check_test_data(void) { + FFTSettings fs; + poly expected_eval, expected_poly, tmp_poly; + new_poly(&expected_eval, 16); + new_poly(&expected_poly, 16); + new_poly(&tmp_poly, 16); + new_fft_settings(&fs, 4); + + for (int i = 0; i < 16; i++) { + fr_from_uint64s(&expected_eval.coeffs[i], expected_eval_u64[i]); + fr_from_uint64s(&expected_poly.coeffs[i], expected_poly_u64[i]); + } + + // Polynomial evalutes to zero at the expected places + for (int i = 0; i < 16; i++) { + if (!exists[i]) { + fr_t tmp; + eval_poly(&tmp, &expected_poly, &fs.expanded_roots_of_unity[i]); + TEST_CHECK(fr_is_zero(&tmp)); + TEST_MSG("Failed for i = %d", i); + } + } + + // This is a curiosity + for (int i = 1; i < 8; i++) { + fr_t tmp; + eval_poly(&tmp, &expected_eval, &fs.expanded_roots_of_unity[i]); + TEST_CHECK(fr_is_zero(&tmp)); + TEST_MSG("Failed for i = %d", i); + } + + // The eval poly is the FFT of the zero poly + TEST_CHECK(C_KZG_OK == fft_fr(tmp_poly.coeffs, expected_eval.coeffs, true, tmp_poly.length, &fs)); + for (int i = 0; i < 16; i++) { + TEST_CHECK(fr_equal(&tmp_poly.coeffs[i], &expected_poly.coeffs[i])); + TEST_MSG("Failed for i = %d", i); + } + + free_poly(&expected_poly); + free_poly(&expected_eval); + free_poly(&tmp_poly); + free_fft_settings(&fs); +} + +void zero_poly_known(void) { + FFTSettings fs; + poly expected_eval, expected_poly, zero_eval, zero_poly; + uint64_t missing[16]; + uint64_t len_missing = 0; + new_poly(&expected_eval, 16); + new_poly(&expected_poly, 16); + new_poly(&zero_eval, 16); + new_poly(&zero_poly, 16); + new_fft_settings(&fs, 4); + + for (int i = 0; i < 16; i++) { + fr_from_uint64s(&expected_eval.coeffs[i], expected_eval_u64[i]); + fr_from_uint64s(&expected_poly.coeffs[i], expected_poly_u64[i]); + if (!exists[i]) { + missing[len_missing++] = i; + } + } + + TEST_CHECK(C_KZG_OK == zero_polynomial_via_multiplication(zero_eval.coeffs, &zero_poly, zero_eval.length, missing, + len_missing, &fs)); + + TEST_CHECK(len_missing + 1 == zero_poly.length); + TEST_MSG("Expected %lu, got %lu", len_missing + 1, zero_poly.length); + + for (int i = 0; i < expected_eval.length; i++) { + TEST_CHECK(fr_equal(&expected_eval.coeffs[i], &zero_eval.coeffs[i])); + TEST_CHECK(fr_equal(&expected_poly.coeffs[i], &zero_poly.coeffs[i])); + } + + free_poly(&expected_poly); + free_poly(&expected_eval); + free_poly(&zero_poly); + free_poly(&zero_eval); + free_fft_settings(&fs); +} + +void zero_poly_random(void) { + for (int its = 0; its < 8; its++) { + srand(its); + for (int scale = 3; scale < 13; scale++) { + FFTSettings fs; + TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, scale)); + + uint64_t *missing; + TEST_CHECK(C_KZG_OK == new_uint64_array(&missing, fs.max_width)); + int len_missing = 0; + + for (int i = 0; i < fs.max_width; i++) { + if (rand() % 2) { + missing[len_missing++] = i; + } + } + + // We know it doesn't work when all indices are missing + if (len_missing == fs.max_width) { + free_fft_settings(&fs); + continue; + } + + fr_t *zero_eval; + poly zero_poly; + TEST_CHECK(C_KZG_OK == new_fr_array(&zero_eval, fs.max_width)); + TEST_CHECK(C_KZG_OK == new_poly(&zero_poly, fs.max_width)); + TEST_CHECK(C_KZG_OK == zero_polynomial_via_multiplication(zero_eval, &zero_poly, fs.max_width, missing, + len_missing, &fs)); + + TEST_CHECK(len_missing + 1 == zero_poly.length); + TEST_MSG("ZeroPolyLen: expected %d, got %lu", len_missing + 1, zero_poly.length); + + int ret = 0; + for (int i = 0; i < len_missing; i++) { + fr_t out; + eval_poly(&out, &zero_poly, &fs.expanded_roots_of_unity[missing[i]]); + ret = TEST_CHECK(fr_is_zero(&out)); + TEST_MSG("Failed for missing[%d] = %lu", i, missing[i]); + } + + fr_t *zero_eval_fft; + TEST_CHECK(C_KZG_OK == new_fr_array(&zero_eval_fft, fs.max_width)); + TEST_CHECK(C_KZG_OK == fft_fr(zero_eval_fft, zero_eval, true, fs.max_width, &fs)); + for (uint64_t i = 0; i < zero_poly.length; i++) { + TEST_CHECK(fr_equal(&zero_poly.coeffs[i], &zero_eval_fft[i])); + } + for (uint64_t i = zero_poly.length; i < fs.max_width; i++) { + TEST_CHECK(fr_is_zero(&zero_eval_fft[i])); + } + + free(missing); + free_poly(&zero_poly); + free(zero_eval); + free(zero_eval_fft); + free_fft_settings(&fs); + } + } +} + +// This didn't work in the original version (ported from Go), but ought now be fine +void zero_poly_all_but_one(void) { + FFTSettings fs; + TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 8)); + + uint64_t *missing; + TEST_CHECK(C_KZG_OK == new_uint64_array(&missing, fs.max_width)); + + // All but the first are missing + for (int i = 0; i < fs.max_width - 1; i++) { + missing[i] = i + 1; + } + int len_missing = fs.max_width - 1; + + fr_t *zero_eval; + poly zero_poly; + TEST_CHECK(C_KZG_OK == new_fr_array(&zero_eval, fs.max_width)); + TEST_CHECK(C_KZG_OK == new_poly(&zero_poly, fs.max_width)); + TEST_CHECK(C_KZG_OK == + zero_polynomial_via_multiplication(zero_eval, &zero_poly, fs.max_width, missing, len_missing, &fs)); + + TEST_CHECK(len_missing + 1 == zero_poly.length); + TEST_MSG("ZeroPolyLen: expected %d, got %lu", len_missing + 1, zero_poly.length); + + int ret = 0; + for (int i = 0; i < len_missing; i++) { + fr_t out; + eval_poly(&out, &zero_poly, &fs.expanded_roots_of_unity[missing[i]]); + ret = TEST_CHECK(fr_is_zero(&out)); + TEST_MSG("Failed for missing[%d] = %lu", i, missing[i]); + } + + fr_t *zero_eval_fft; + TEST_CHECK(C_KZG_OK == new_fr_array(&zero_eval_fft, fs.max_width)); + TEST_CHECK(C_KZG_OK == fft_fr(zero_eval_fft, zero_eval, true, fs.max_width, &fs)); + for (uint64_t i = 0; i < zero_poly.length; i++) { + TEST_CHECK(fr_equal(&zero_poly.coeffs[i], &zero_eval_fft[i])); + } + for (uint64_t i = zero_poly.length; i < fs.max_width; i++) { + TEST_CHECK(fr_is_zero(&zero_eval_fft[i])); + } + + free(missing); + free_poly(&zero_poly); + free(zero_eval); + free(zero_eval_fft); + free_fft_settings(&fs); +} + +// This is to prevent regressions - 252 missing at width 8 is an edge case which has 4 full partials +void zero_poly_252(void) { + FFTSettings fs; + TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 8)); + + uint64_t *missing; + TEST_CHECK(C_KZG_OK == new_uint64_array(&missing, fs.max_width)); + + int len_missing = 252; + for (int i = 0; i < len_missing; i++) { + missing[i] = i; + } + + fr_t *zero_eval; + poly zero_poly; + TEST_CHECK(C_KZG_OK == new_fr_array(&zero_eval, fs.max_width)); + TEST_CHECK(C_KZG_OK == new_poly(&zero_poly, fs.max_width)); + TEST_CHECK(C_KZG_OK == + zero_polynomial_via_multiplication(zero_eval, &zero_poly, fs.max_width, missing, len_missing, &fs)); + + TEST_CHECK(len_missing + 1 == zero_poly.length); + TEST_MSG("ZeroPolyLen: expected %d, got %lu", len_missing + 1, zero_poly.length); + + int ret = 0; + for (int i = 0; i < len_missing; i++) { + fr_t out; + eval_poly(&out, &zero_poly, &fs.expanded_roots_of_unity[missing[i]]); + ret = TEST_CHECK(fr_is_zero(&out)); + TEST_MSG("Failed for missing[%d] = %lu", i, missing[i]); + } + + fr_t *zero_eval_fft; + TEST_CHECK(C_KZG_OK == new_fr_array(&zero_eval_fft, fs.max_width)); + TEST_CHECK(C_KZG_OK == fft_fr(zero_eval_fft, zero_eval, true, fs.max_width, &fs)); + for (uint64_t i = 0; i < zero_poly.length; i++) { + TEST_CHECK(fr_equal(&zero_poly.coeffs[i], &zero_eval_fft[i])); + } + for (uint64_t i = zero_poly.length; i < fs.max_width; i++) { + TEST_CHECK(fr_is_zero(&zero_eval_fft[i])); + } + + free(missing); + free_poly(&zero_poly); + free(zero_eval); + free(zero_eval_fft); + free_fft_settings(&fs); +} + +TEST_LIST = { + {"ZERO_POLY_TEST", title}, + {"test_reduce_partials", test_reduce_partials}, + {"check_test_data", check_test_data}, + {"reduce_partials_random", reduce_partials_random}, + {"zero_poly_known", zero_poly_known}, + {"zero_poly_random", zero_poly_random}, + {"zero_poly_all_but_one", zero_poly_all_but_one}, + {"zero_poly_252", zero_poly_252}, + {NULL, NULL} /* zero record marks the end of the list */ +}; + +#endif // KZGTEST diff --git a/src/zero_poly.h b/src/zero_poly.h index 64d0358..0c78567 100644 --- a/src/zero_poly.h +++ b/src/zero_poly.h @@ -24,10 +24,6 @@ #include "fft_common.h" #include "poly.h" -C_KZG_RET do_zero_poly_mul_partial(poly *dst, const uint64_t *indices, uint64_t len_indices, uint64_t stride, - const FFTSettings *fs); -C_KZG_RET reduce_partials(poly *dst, uint64_t len_dst, fr_t *scratch, uint64_t len_scratch, const poly *partials, - uint64_t partial_count, const FFTSettings *fs); C_KZG_RET zero_polynomial_via_multiplication(fr_t *zero_eval, poly *zero_poly, uint64_t width, const uint64_t *missing_indices, uint64_t len_missing, const FFTSettings *fs); diff --git a/src/zero_poly_test.c b/src/zero_poly_test.c deleted file mode 100644 index 12834bc..0000000 --- a/src/zero_poly_test.c +++ /dev/null @@ -1,414 +0,0 @@ -/* - * Copyright 2021 Benjamin Edgington - * - * Licensed under the Apache License, Version 2.0 (the "License"); - * you may not use this file except in compliance with the License. - * You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -#include "../inc/acutest.h" -#include "c_kzg_util.h" -#include "test_util.h" -#include "zero_poly.h" -#include "fft_fr.h" -#include "debug_util.h" - -bool exists[16] = {true, false, false, true, false, true, true, false, - false, false, true, true, false, true, false, true}; - -uint64_t expected_eval_u64[16][4] = { - {0xfd5a5130b97ce0c3L, 0xb4748a4cb0f90e6dL, 0x12a1ab34b25b18c1L, 0x5a5ac0c81c9f7ea8L}, - {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, - {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, - {0xaa385cbce3dd1657L, 0x2fdab57a38bdb514L, 0x20e022e205dafa53L, 0x14077dd3f5d996b1L}, - {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, - {0x194018614b6f7276L, 0xdf2b18f870532376L, 0x1ff427cd5b583fe6L, 0x014d6444ff03dd09L}, - {0xcc84c2de684c0ddeL, 0xf1e7ab32aa830d02L, 0x967bf35a2a691f20L, 0x046109731cdf0d3cL}, - {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, - {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, - {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, - {0x96cddd2924212afbL, 0xeaa4c1f51421d8d8L, 0x3ae969cfa34d0ed1L, 0x6b6c5e876bc3916dL}, - {0x449310802f74ad49L, 0x47c940979163037aL, 0x10d311564afb9b2aL, 0x269b8531c369bafbL}, - {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, - {0xd9af75fe35c16cf1L, 0x068bb140cea92f75L, 0xe769811965e10a47L, 0x48ed97e6745612f2L}, - {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, - {0x7ef1f59bb1677685L, 0x33a637296680e8ceL, 0xaaf62b3f6e016709L, 0x454a299178a4dba9L}}; - -uint64_t expected_poly_u64[16][4] = { - {0xac159e2688bd4333L, 0x3bfef0f00df2ec88L, 0x561dcd0fd4d314d9L, 0x533bd8c1e977024eL}, - {0x18bc6eedc010ef8dL, 0xc731a3eb4ea2ab70L, 0x5c2589357ae121a8L, 0x04f9108d308f7016L}, - {0x232759f49556ac08L, 0x9776fe2e9f4c613cL, 0x74d5bed4eb2de960L, 0x1f6cf6719bfa0e68L}, - {0xf2f3461e8ab1ae34L, 0xeb220fcc11ef1c80L, 0x7a4637d3a637739bL, 0x19901a58cd177c53L}, - {0x9340f62465a1f4feL, 0xd9cb3ea6de494a11L, 0xee92ebc763cdff5dL, 0x5443e89811b5b9f5L}, - {0x269a255e2e4e48a4L, 0xfadae7a89d9b2f2bL, 0xb5515799b41e1a88L, 0x2e990979a0ffcee5L}, - {0x1c2f3a5759088c29L, 0x2a958d654cf1795fL, 0x9ca121fa43d152d1L, 0x1425239535953093L}, - {0x4c634e2d63ad89fdL, 0xd6ea7bc7da4ebe1aL, 0x9730a8fb88c7c895L, 0x1a01ffae0477c2a8L}, - {0x0000000000000001L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, - {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, - {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, - {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, - {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, - {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, - {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}, - {0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L}}; - -void test_reduce_partials(void) { - FFTSettings fs; - TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 4)); - fr_t from_tree_reduction_coeffs[16], from_direct_coeffs[9], scratch[48]; - poly from_tree_reduction, from_direct; - from_tree_reduction.coeffs = from_tree_reduction_coeffs; - from_tree_reduction.length = 16; - from_direct.coeffs = from_direct_coeffs; - from_direct.length = 9; - - // Via reduce_partials - - poly partials[4]; - fr_t partial0[3], partial1[3], partial2[3], partial3[3]; - partials[0].coeffs = partial0, partials[0].length = 3; - partials[1].coeffs = partial1, partials[1].length = 3; - partials[2].coeffs = partial2, partials[2].length = 3; - partials[3].coeffs = partial3, partials[3].length = 3; - const uint64_t partial_indices[4][2] = {{1, 3}, {7, 8}, {9, 10}, {12, 13}}; - for (int i = 0; i < 4; i++) { - TEST_CHECK(C_KZG_OK == do_zero_poly_mul_partial(&partials[i], partial_indices[i], 2, 1, &fs)); - } - TEST_CHECK(C_KZG_OK == reduce_partials(&from_tree_reduction, 16, scratch, 48, partials, 4, &fs)); - - // Direct - uint64_t indices[] = {1, 3, 7, 8, 9, 10, 12, 13}; - TEST_CHECK(C_KZG_OK == do_zero_poly_mul_partial(&from_direct, indices, 8, 1, &fs)); - - // Compare - for (int i = 0; i < 9; i++) { - TEST_CHECK(fr_equal(&from_tree_reduction.coeffs[i], &from_direct.coeffs[i])); - } - - free_fft_settings(&fs); -} - -void reduce_partials_random(void) { - for (int scale = 5; scale < 13; scale++) { - for (int ii = 1; ii <= 7; ii++) { - float missing_ratio = 0.1 * ii; - - FFTSettings fs; - new_fft_settings(&fs, scale); - uint64_t point_count = fs.max_width; - uint64_t missing_count = point_count * missing_ratio; - - uint64_t *missing; - TEST_CHECK(C_KZG_OK == new_uint64_array(&missing, point_count)); - for (uint64_t i = 0; i < point_count; i++) { - missing[i] = i; - } - shuffle(missing, point_count); - - // Build the partials - poly *partials; - const int missing_per_partial = 63; - uint64_t indices[missing_per_partial]; - uint64_t partial_count = (missing_count + missing_per_partial - 1) / missing_per_partial; - TEST_CHECK(C_KZG_OK == new_poly_array(&partials, partial_count)); - for (uint64_t i = 0; i < partial_count; i++) { - uint64_t start = i * missing_per_partial; - uint64_t end = start + missing_per_partial; - if (end > missing_count) end = missing_count; - uint64_t partial_size = end - start; - TEST_CHECK(C_KZG_OK == new_fr_array(&partials[i].coeffs, partial_size + 1)); - for (int j = 0; j < partial_size; j++) { - indices[j] = missing[i * missing_per_partial + j]; - } - partials[i].length = partial_size + 1; - TEST_CHECK(C_KZG_OK == do_zero_poly_mul_partial(&partials[i], indices, partial_size, 1, &fs)); - } - - // From tree reduction - poly from_tree_reduction; - TEST_CHECK(C_KZG_OK == new_poly(&from_tree_reduction, point_count)); - fr_t *scratch; - TEST_CHECK(C_KZG_OK == new_fr_array(&scratch, point_count * 3)); - TEST_CHECK(C_KZG_OK == reduce_partials(&from_tree_reduction, point_count, scratch, point_count * 3, - partials, partial_count, &fs)); - - // From direct - poly from_direct; - TEST_CHECK(C_KZG_OK == new_poly(&from_direct, missing_count + 1)); - TEST_CHECK(C_KZG_OK == - do_zero_poly_mul_partial(&from_direct, missing, missing_count, fs.max_width / point_count, &fs)); - - for (uint64_t i = 0; i < missing_count + 1; i++) { - TEST_CHECK(fr_equal(&from_tree_reduction.coeffs[i], &from_direct.coeffs[i])); - } - - free_poly(&from_tree_reduction); - free_poly(&from_direct); - free(scratch); - for (uint64_t i = 0; i < partial_count; i++) { - free_poly(&partials[i]); - } - free(partials); - free(missing); - free_fft_settings(&fs); - } - } -} - -void check_test_data(void) { - FFTSettings fs; - poly expected_eval, expected_poly, tmp_poly; - new_poly(&expected_eval, 16); - new_poly(&expected_poly, 16); - new_poly(&tmp_poly, 16); - new_fft_settings(&fs, 4); - - for (int i = 0; i < 16; i++) { - fr_from_uint64s(&expected_eval.coeffs[i], expected_eval_u64[i]); - fr_from_uint64s(&expected_poly.coeffs[i], expected_poly_u64[i]); - } - - // Polynomial evalutes to zero at the expected places - for (int i = 0; i < 16; i++) { - if (!exists[i]) { - fr_t tmp; - eval_poly(&tmp, &expected_poly, &fs.expanded_roots_of_unity[i]); - TEST_CHECK(fr_is_zero(&tmp)); - TEST_MSG("Failed for i = %d", i); - } - } - - // This is a curiosity - for (int i = 1; i < 8; i++) { - fr_t tmp; - eval_poly(&tmp, &expected_eval, &fs.expanded_roots_of_unity[i]); - TEST_CHECK(fr_is_zero(&tmp)); - TEST_MSG("Failed for i = %d", i); - } - - // The eval poly is the FFT of the zero poly - TEST_CHECK(C_KZG_OK == fft_fr(tmp_poly.coeffs, expected_eval.coeffs, true, tmp_poly.length, &fs)); - for (int i = 0; i < 16; i++) { - TEST_CHECK(fr_equal(&tmp_poly.coeffs[i], &expected_poly.coeffs[i])); - TEST_MSG("Failed for i = %d", i); - } - - free_poly(&expected_poly); - free_poly(&expected_eval); - free_poly(&tmp_poly); - free_fft_settings(&fs); -} - -void zero_poly_known(void) { - FFTSettings fs; - poly expected_eval, expected_poly, zero_eval, zero_poly; - uint64_t missing[16]; - uint64_t len_missing = 0; - new_poly(&expected_eval, 16); - new_poly(&expected_poly, 16); - new_poly(&zero_eval, 16); - new_poly(&zero_poly, 16); - new_fft_settings(&fs, 4); - - for (int i = 0; i < 16; i++) { - fr_from_uint64s(&expected_eval.coeffs[i], expected_eval_u64[i]); - fr_from_uint64s(&expected_poly.coeffs[i], expected_poly_u64[i]); - if (!exists[i]) { - missing[len_missing++] = i; - } - } - - TEST_CHECK(C_KZG_OK == zero_polynomial_via_multiplication(zero_eval.coeffs, &zero_poly, zero_eval.length, missing, - len_missing, &fs)); - - TEST_CHECK(len_missing + 1 == zero_poly.length); - TEST_MSG("Expected %lu, got %lu", len_missing + 1, zero_poly.length); - - for (int i = 0; i < expected_eval.length; i++) { - TEST_CHECK(fr_equal(&expected_eval.coeffs[i], &zero_eval.coeffs[i])); - TEST_CHECK(fr_equal(&expected_poly.coeffs[i], &zero_poly.coeffs[i])); - } - - free_poly(&expected_poly); - free_poly(&expected_eval); - free_poly(&zero_poly); - free_poly(&zero_eval); - free_fft_settings(&fs); -} - -void zero_poly_random(void) { - for (int its = 0; its < 8; its++) { - srand(its); - for (int scale = 3; scale < 13; scale++) { - FFTSettings fs; - TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, scale)); - - uint64_t *missing; - TEST_CHECK(C_KZG_OK == new_uint64_array(&missing, fs.max_width)); - int len_missing = 0; - - for (int i = 0; i < fs.max_width; i++) { - if (rand() % 2) { - missing[len_missing++] = i; - } - } - - // We know it doesn't work when all indices are missing - if (len_missing == fs.max_width) { - free_fft_settings(&fs); - continue; - } - - fr_t *zero_eval; - poly zero_poly; - TEST_CHECK(C_KZG_OK == new_fr_array(&zero_eval, fs.max_width)); - TEST_CHECK(C_KZG_OK == new_poly(&zero_poly, fs.max_width)); - TEST_CHECK(C_KZG_OK == zero_polynomial_via_multiplication(zero_eval, &zero_poly, fs.max_width, missing, - len_missing, &fs)); - - TEST_CHECK(len_missing + 1 == zero_poly.length); - TEST_MSG("ZeroPolyLen: expected %d, got %lu", len_missing + 1, zero_poly.length); - - int ret = 0; - for (int i = 0; i < len_missing; i++) { - fr_t out; - eval_poly(&out, &zero_poly, &fs.expanded_roots_of_unity[missing[i]]); - ret = TEST_CHECK(fr_is_zero(&out)); - TEST_MSG("Failed for missing[%d] = %lu", i, missing[i]); - } - - fr_t *zero_eval_fft; - TEST_CHECK(C_KZG_OK == new_fr_array(&zero_eval_fft, fs.max_width)); - TEST_CHECK(C_KZG_OK == fft_fr(zero_eval_fft, zero_eval, true, fs.max_width, &fs)); - for (uint64_t i = 0; i < zero_poly.length; i++) { - TEST_CHECK(fr_equal(&zero_poly.coeffs[i], &zero_eval_fft[i])); - } - for (uint64_t i = zero_poly.length; i < fs.max_width; i++) { - TEST_CHECK(fr_is_zero(&zero_eval_fft[i])); - } - - free(missing); - free_poly(&zero_poly); - free(zero_eval); - free(zero_eval_fft); - free_fft_settings(&fs); - } - } -} - -// This didn't work in the original version (ported from Go), but ought now be fine -void zero_poly_all_but_one(void) { - FFTSettings fs; - TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 8)); - - uint64_t *missing; - TEST_CHECK(C_KZG_OK == new_uint64_array(&missing, fs.max_width)); - - // All but the first are missing - for (int i = 0; i < fs.max_width - 1; i++) { - missing[i] = i + 1; - } - int len_missing = fs.max_width - 1; - - fr_t *zero_eval; - poly zero_poly; - TEST_CHECK(C_KZG_OK == new_fr_array(&zero_eval, fs.max_width)); - TEST_CHECK(C_KZG_OK == new_poly(&zero_poly, fs.max_width)); - TEST_CHECK(C_KZG_OK == - zero_polynomial_via_multiplication(zero_eval, &zero_poly, fs.max_width, missing, len_missing, &fs)); - - TEST_CHECK(len_missing + 1 == zero_poly.length); - TEST_MSG("ZeroPolyLen: expected %d, got %lu", len_missing + 1, zero_poly.length); - - int ret = 0; - for (int i = 0; i < len_missing; i++) { - fr_t out; - eval_poly(&out, &zero_poly, &fs.expanded_roots_of_unity[missing[i]]); - ret = TEST_CHECK(fr_is_zero(&out)); - TEST_MSG("Failed for missing[%d] = %lu", i, missing[i]); - } - - fr_t *zero_eval_fft; - TEST_CHECK(C_KZG_OK == new_fr_array(&zero_eval_fft, fs.max_width)); - TEST_CHECK(C_KZG_OK == fft_fr(zero_eval_fft, zero_eval, true, fs.max_width, &fs)); - for (uint64_t i = 0; i < zero_poly.length; i++) { - TEST_CHECK(fr_equal(&zero_poly.coeffs[i], &zero_eval_fft[i])); - } - for (uint64_t i = zero_poly.length; i < fs.max_width; i++) { - TEST_CHECK(fr_is_zero(&zero_eval_fft[i])); - } - - free(missing); - free_poly(&zero_poly); - free(zero_eval); - free(zero_eval_fft); - free_fft_settings(&fs); -} - -// This is to prevent regressions - 252 missing at width 8 is an edge case which has 4 full partials -void zero_poly_252(void) { - FFTSettings fs; - TEST_CHECK(C_KZG_OK == new_fft_settings(&fs, 8)); - - uint64_t *missing; - TEST_CHECK(C_KZG_OK == new_uint64_array(&missing, fs.max_width)); - - int len_missing = 252; - for (int i = 0; i < len_missing; i++) { - missing[i] = i; - } - - fr_t *zero_eval; - poly zero_poly; - TEST_CHECK(C_KZG_OK == new_fr_array(&zero_eval, fs.max_width)); - TEST_CHECK(C_KZG_OK == new_poly(&zero_poly, fs.max_width)); - TEST_CHECK(C_KZG_OK == - zero_polynomial_via_multiplication(zero_eval, &zero_poly, fs.max_width, missing, len_missing, &fs)); - - TEST_CHECK(len_missing + 1 == zero_poly.length); - TEST_MSG("ZeroPolyLen: expected %d, got %lu", len_missing + 1, zero_poly.length); - - int ret = 0; - for (int i = 0; i < len_missing; i++) { - fr_t out; - eval_poly(&out, &zero_poly, &fs.expanded_roots_of_unity[missing[i]]); - ret = TEST_CHECK(fr_is_zero(&out)); - TEST_MSG("Failed for missing[%d] = %lu", i, missing[i]); - } - - fr_t *zero_eval_fft; - TEST_CHECK(C_KZG_OK == new_fr_array(&zero_eval_fft, fs.max_width)); - TEST_CHECK(C_KZG_OK == fft_fr(zero_eval_fft, zero_eval, true, fs.max_width, &fs)); - for (uint64_t i = 0; i < zero_poly.length; i++) { - TEST_CHECK(fr_equal(&zero_poly.coeffs[i], &zero_eval_fft[i])); - } - for (uint64_t i = zero_poly.length; i < fs.max_width; i++) { - TEST_CHECK(fr_is_zero(&zero_eval_fft[i])); - } - - free(missing); - free_poly(&zero_poly); - free(zero_eval); - free(zero_eval_fft); - free_fft_settings(&fs); -} - -TEST_LIST = { - {"ZERO_POLY_TEST", title}, - {"test_reduce_partials", test_reduce_partials}, - {"check_test_data", check_test_data}, - {"reduce_partials_random", reduce_partials_random}, - {"zero_poly_known", zero_poly_known}, - {"zero_poly_random", zero_poly_random}, - {"zero_poly_all_but_one", zero_poly_all_but_one}, - {"zero_poly_252", zero_poly_252}, - {NULL, NULL} /* zero record marks the end of the list */ -};