mirror of
https://github.com/status-im/c-kzg-4844.git
synced 2025-02-24 07:48:05 +00:00
Split out FFT utilities and tests to separate files
This commit is contained in:
parent
8de064215a
commit
3ee603fbdc
@ -1,10 +1,14 @@
|
||||
tests = fft_fr_test
|
||||
tests = fft_fr_test fft_util_test
|
||||
|
||||
.PRECIOUS: %.o
|
||||
|
||||
%.o: %.c %.h c-kzg.h
|
||||
clang -Wall -c $*.c
|
||||
|
||||
fft_fr_test: fft_fr.o fft_fr_test.c test_util.o fft_util.o
|
||||
clang -Wall -o $@ $@.c test_util.o fft_fr.o fft_util.o -L../lib -lblst
|
||||
./$@
|
||||
|
||||
%_test: %.o %_test.c test_util.o
|
||||
clang -Wall -o $@ $@.c test_util.o $*.o -L../lib -lblst
|
||||
./$@
|
||||
|
@ -15,4 +15,5 @@
|
||||
*/
|
||||
|
||||
#include <stdbool.h>
|
||||
#include <assert.h>
|
||||
#include "../inc/blst.h"
|
||||
|
56
src/fft_fr.c
56
src/fft_fr.c
@ -16,62 +16,6 @@
|
||||
|
||||
#include "fft_fr.h"
|
||||
|
||||
bool is_one(const blst_fr *fr_p) {
|
||||
uint64_t a[4];
|
||||
blst_uint64_from_fr(a, fr_p);
|
||||
return a[0] == 1 && a[1] == 0 && a[2] == 0 && a[3] == 0;
|
||||
}
|
||||
|
||||
bool is_power_of_two(uint64_t n) {
|
||||
return (n & (n - 1)) == 0;
|
||||
}
|
||||
|
||||
void fr_from_uint64(blst_fr *a, uint64_t n) {
|
||||
uint64_t vals[] = {n, 0, 0, 0};
|
||||
blst_fr_from_uint64(a, vals);
|
||||
}
|
||||
|
||||
// Returns an array of powers of the root of unity
|
||||
// Allocates space for the array that needs to be freed later
|
||||
blst_fr *expand_root_of_unity(blst_fr *root_of_unity, uint64_t width) {
|
||||
blst_fr *roots = malloc((width + 1) * sizeof(blst_fr));
|
||||
roots[0] = one;
|
||||
roots[1] = *root_of_unity;
|
||||
|
||||
for (int i = 2; !is_one(&roots[i - 1]); i++) {
|
||||
assert(i <= width);
|
||||
blst_fr_mul(&roots[i], &roots[i - 1], root_of_unity);
|
||||
}
|
||||
assert(is_one(&roots[width]));
|
||||
|
||||
return roots;
|
||||
}
|
||||
|
||||
// Return a reversed copy of the list of Fr provided
|
||||
// `width` is one less than the length of `r`
|
||||
// Allocates space for the array that needs to be freed later
|
||||
blst_fr *reverse(blst_fr *r, uint64_t width) {
|
||||
blst_fr *rr = malloc((width + 1) * sizeof(blst_fr));
|
||||
for (int i = 0; i <= width; i++) {
|
||||
rr[i] = r[width - i];
|
||||
}
|
||||
return rr;
|
||||
}
|
||||
|
||||
FFTSettings new_fft_settings(unsigned int max_scale) {
|
||||
FFTSettings s;
|
||||
s.max_width = (uint64_t)1 << max_scale;
|
||||
blst_fr_from_uint64(&s.root_of_unity, scale2_root_of_unity[max_scale]);
|
||||
s.expanded_roots_of_unity = expand_root_of_unity(&s.root_of_unity, s.max_width);
|
||||
s.reverse_roots_of_unity = reverse(s.expanded_roots_of_unity, s.max_width);
|
||||
return s;
|
||||
}
|
||||
|
||||
void free_fft_settings(FFTSettings *s) {
|
||||
free(s->expanded_roots_of_unity);
|
||||
free(s->reverse_roots_of_unity);
|
||||
}
|
||||
|
||||
// Slow Fourier Transform (simple, good for small sizes)
|
||||
void slow_ft(blst_fr *out, blst_fr *in, uint64_t offset, uint64_t stride, blst_fr *roots, uint64_t roots_stride, uint64_t l) {
|
||||
blst_fr v, last, tmp;
|
||||
|
63
src/fft_fr.h
63
src/fft_fr.h
@ -14,70 +14,9 @@
|
||||
* limitations under the License.
|
||||
*/
|
||||
|
||||
#include <stdlib.h>
|
||||
#include <assert.h>
|
||||
#include <stdio.h>
|
||||
#include "c-kzg.h"
|
||||
#include "fft_util.h"
|
||||
|
||||
// This is 1 in Blst's `blst_fr` limb representation. Crazy but true.
|
||||
static const blst_fr one =
|
||||
{0x00000001fffffffeL, 0x5884b7fa00034802L, 0x998c4fefecbc4ff5L, 0x1824b159acc5056fL};
|
||||
|
||||
// MODULUS = 52435875175126190479447740508185965837690552500527637822603658699938581184513
|
||||
// PRIMITIVE_ROOT = 5
|
||||
// [pow(PRIMITIVE_ROOT, (MODULUS - 1) // (2**i), MODULUS) for i in range(32)]
|
||||
//
|
||||
// These are not in `blst_fr` limb format and must be converted via `blst_fr_from_uint64()`
|
||||
static const uint64_t scale2_root_of_unity[][4] =
|
||||
{
|
||||
{0x0000000000000001L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L},
|
||||
{0xffffffff00000000L, 0x53bda402fffe5bfeL, 0x3339d80809a1d805L, 0x73eda753299d7d48L},
|
||||
{0x0001000000000000L, 0xec03000276030000L, 0x8d51ccce760304d0L, 0x0000000000000000L},
|
||||
{0x8dd702cb688bc087L, 0xa032824078eaa4feL, 0xa733b23a98ca5b22L, 0x3f96405d25a31660L},
|
||||
{0x95df4b360411fe73L, 0x1ef4e672ebc1e1bbL, 0x6e92a9c460afca4aL, 0x4f2c596e753e4fccL},
|
||||
{0x9733f3a6ba60eaa6L, 0xbd7fdf9c77487ae7L, 0xd84f8612c8b6cc00L, 0x476fa2fb6162ffabL},
|
||||
{0xd2fc5e69ac5db47fL, 0xa12a70a615d0b8e4L, 0x293b1d67bc8de5d9L, 0x0e4840ac57f86f5eL},
|
||||
{0xb750da4cab28e208L, 0x501dff643be95635L, 0x8cbe2437f0b4b276L, 0x07d0c802a94a946eL},
|
||||
{0x2cabadec2fe322b8L, 0x752c84f315412560L, 0x32a732ae1a3b0aefL, 0x2e95da59a33dcbf2L},
|
||||
{0x33811ea1fe0c65f4L, 0x15c1ad4c687f28a2L, 0xecfbede342dee7f4L, 0x1bb466679a5d88b1L},
|
||||
{0xd58a5af42d010ff9L, 0x79efd6b0570bf109L, 0x3ed6d55a6350721dL, 0x2f27b09858f43cefL},
|
||||
{0x74a1f6718c130477L, 0xa534af14b61e0abeL, 0xeb674a1a620890d7L, 0x43527a8bca252472L},
|
||||
{0x450d9f977ea8ee05L, 0x565af17137d56fc0L, 0xe155cb4893f9e9acL, 0x110cebd0c8e9101bL},
|
||||
{0x23c9159959a0be92L, 0x87d188ce7a027759L, 0x70491431cab3c3ccL, 0x0ac00eb8b3f7f8daL},
|
||||
{0x13e96ade69583404L, 0x82c057275306243dL, 0x77e48bf529ca9f2aL, 0x50646ac81fe19595L},
|
||||
{0xe6a354dda97eccd4L, 0x39929d2e88fbbc57L, 0xa22ba63dd6e7b1c8L, 0x42c22911f5f07f43L},
|
||||
{0x137b458acfc35f7aL, 0x0caba63a29c01b06L, 0x0409ee987a02402cL, 0x6709c6cd56aa725bL},
|
||||
{0x10251f7d8831e03eL, 0x77d85a937ff858ecL, 0xebe905bd4fb9ac5cL, 0x05deb333f8727901L},
|
||||
{0xbf87b689b9009408L, 0x4f730e7ddd3ccc96L, 0xfd7f05ba4610300cL, 0x5ef5e8db0b8ac903L},
|
||||
{0x6499688417cd0c14L, 0xa672867368812f7fL, 0x2e1d9a1922cc3253L, 0x3a689e83aa0a1d80L},
|
||||
{0x20b53cbe41144deaL, 0x870c46fac2f0fcbdL, 0x556c35f6537d6971L, 0x3436287f5f686d91L},
|
||||
{0x007e082a436ba2e7L, 0x67c6630f9116e877L, 0x36f8f165fb4460f7L, 0x6eee34d57e7046e0L},
|
||||
{0xc5b670eea53a56d1L, 0x127d1f4253037d7bL, 0x57d4257ea722c2e2L, 0x03ae26a333cbd838L},
|
||||
{0x1e91484876504cf8L, 0x55bbbf1eb63edd02L, 0xbcdafec84e55aa02L, 0x5145c4cd2dc0beb0L},
|
||||
{0x5b90153a1ab70e2cL, 0x8deffa3175fb0ab8L, 0xc553ae2346900c95L, 0x1d31dcdc6bd3118cL},
|
||||
{0x801c894c59a2e8ebL, 0xbc535c5ce12fc974L, 0x95508d2747d39803L, 0x16d9d3cdac5d094fL},
|
||||
{0x810fa372cca1d8beL, 0xc67b8c2882e0bfa7L, 0xdbb4edf0e2d35bc2L, 0x712d15805087c995L},
|
||||
{0xeb162203fd88f133L, 0xac96c38ff010ea74L, 0x4307987fe64cfc70L, 0x350fe98d37b7a114L},
|
||||
{0xaba2f51842f2a254L, 0x4d7f3c3aa71efc0cL, 0x97ae418dd274a80aL, 0x2967385d5e3e7682L},
|
||||
{0x75c55c7b575a0b79L, 0x3ba4a15774a7ded1L, 0xc3974d73a04fccf3L, 0x705aba4f4a939684L},
|
||||
{0x8409a9ea14ebb608L, 0xfad0084e66bac611L, 0x04287254811c1dfbL, 0x086d072b23b30c29L},
|
||||
{0xb427c9b367e4756aL, 0xc7537fb902ebc38dL, 0x51de21becd6a205fL, 0x6064ab727923597dL}
|
||||
};
|
||||
|
||||
typedef struct {
|
||||
uint64_t max_width;
|
||||
blst_fr root_of_unity;
|
||||
blst_fr *expanded_roots_of_unity;
|
||||
blst_fr *reverse_roots_of_unity;
|
||||
} FFTSettings;
|
||||
|
||||
bool is_one(const blst_fr *fr_p);
|
||||
bool is_power_of_two(uint64_t n);
|
||||
void fr_from_uint64(blst_fr *a, uint64_t n);
|
||||
blst_fr *expand_root_of_unity(blst_fr *root_of_unity, uint64_t width);
|
||||
blst_fr *reverse(blst_fr *r, uint64_t width);
|
||||
FFTSettings new_fft_settings(unsigned int max_scale);
|
||||
void free_fft_settings(FFTSettings *s);
|
||||
void slow_ft(blst_fr *out, blst_fr *in, uint64_t offset, uint64_t stride, blst_fr *roots, uint64_t roots_stride, uint64_t l);
|
||||
void fast_ft(blst_fr *out, blst_fr *in, uint64_t offset, uint64_t stride, blst_fr *roots, uint64_t roots_stride, uint64_t l);
|
||||
void fft_helper(blst_fr *out, blst_fr *in, uint64_t offset, uint64_t stride, blst_fr *roots, uint64_t roots_stride, uint64_t l);
|
||||
|
@ -18,8 +18,6 @@
|
||||
#include "test_util.h"
|
||||
#include "fft_fr.h"
|
||||
|
||||
#define NUM_ROOTS 32
|
||||
|
||||
const uint64_t inv_fft_expected[][4] =
|
||||
{
|
||||
{0x7fffffff80000008L, 0xa9ded2017fff2dffL, 0x199cec0404d0ec02L, 0x39f6d3a994cebea4L},
|
||||
@ -40,117 +38,6 @@ const uint64_t inv_fft_expected[][4] =
|
||||
{0xd8cda22e753b3117L, 0x880454ec72238f55L, 0xcaf6199fc14a5353L, 0x197df7c2f05866d4L}
|
||||
};
|
||||
|
||||
void is_one_works(void) {
|
||||
TEST_CHECK(true == is_one(&one));
|
||||
}
|
||||
|
||||
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_are_plausible(void) {
|
||||
blst_fr r;
|
||||
for (unsigned int i = 0; i < NUM_ROOTS; i++) {
|
||||
blst_fr_from_uint64(&r, scale2_root_of_unity[i]);
|
||||
for (unsigned int j = 0; j < i; j++) {
|
||||
blst_fr_sqr(&r, &r);
|
||||
}
|
||||
TEST_CHECK(true == is_one(&r));
|
||||
TEST_MSG("Root %d failed", i);
|
||||
}
|
||||
}
|
||||
|
||||
void reverse_works(void) {
|
||||
int n = 24;
|
||||
blst_fr arr[n + 1];
|
||||
blst_fr *rev, diff;
|
||||
|
||||
// Initialise - increasing values
|
||||
arr[0] = one;
|
||||
for (int i = 1; i <= n; i++) {
|
||||
blst_fr_add(arr + i, arr + i - 1, &one);
|
||||
}
|
||||
|
||||
// Reverse
|
||||
rev = reverse(arr, n);
|
||||
|
||||
// Verify - decreasing values
|
||||
for (int i = 0; i < n; i++) {
|
||||
blst_fr_sub(&diff, rev + i, rev + i + 1);
|
||||
TEST_CHECK(true == is_one(&diff));
|
||||
}
|
||||
TEST_CHECK(true == is_one(rev + n));
|
||||
|
||||
free(rev);
|
||||
}
|
||||
|
||||
void expand_roots_is_plausible(void) {
|
||||
// Just test one (largeish) value of scale
|
||||
unsigned int scale = 20;
|
||||
unsigned int width = 1 << scale;
|
||||
blst_fr root, *expanded, prod;
|
||||
|
||||
// Initialise
|
||||
blst_fr_from_uint64(&root, scale2_root_of_unity[scale]);
|
||||
expanded = expand_root_of_unity(&root, width);
|
||||
|
||||
// Verify - each pair should multiply to one
|
||||
TEST_CHECK(true == is_one(expanded + 0));
|
||||
TEST_CHECK(true == is_one(expanded + width));
|
||||
for (unsigned int i = 1; i <= width / 2; i++) {
|
||||
blst_fr_mul(&prod, expanded + i, expanded + width - i);
|
||||
TEST_CHECK(true == is_one(&prod));
|
||||
}
|
||||
|
||||
free(expanded);
|
||||
}
|
||||
|
||||
void new_fft_settings_is_plausible(void) {
|
||||
// Just test one (largeish) value of scale
|
||||
unsigned int scale = 21;
|
||||
unsigned int width = 1 << scale;
|
||||
blst_fr prod;
|
||||
FFTSettings s = new_fft_settings(scale);
|
||||
|
||||
// Verify - each pair should multiply to one
|
||||
for (unsigned int i = 1; i <= width; i++) {
|
||||
blst_fr_mul(&prod, s.expanded_roots_of_unity + i, s.reverse_roots_of_unity + i);
|
||||
TEST_CHECK(true == is_one(&prod));
|
||||
}
|
||||
|
||||
free_fft_settings(&s);
|
||||
}
|
||||
|
||||
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 fr_from_uint64_works(void) {
|
||||
blst_fr a;
|
||||
fr_from_uint64(&a, 1);
|
||||
TEST_CHECK(true == is_one(&a));
|
||||
}
|
||||
|
||||
void fr_equal_works(void) {
|
||||
blst_fr a, b;
|
||||
blst_fr_from_uint64(&a, scale2_root_of_unity[15]);
|
||||
blst_fr_from_uint64(&b, scale2_root_of_unity[16]);
|
||||
TEST_CHECK(true == fr_equal(&a, &a));
|
||||
TEST_CHECK(false == fr_equal(&a, &b));
|
||||
}
|
||||
|
||||
void compare_sft_fft(void) {
|
||||
// Initialise: ascending values of i (could be anything), and arbitrary size
|
||||
unsigned int size = 8;
|
||||
@ -220,15 +107,6 @@ void inverse_fft(void) {
|
||||
|
||||
TEST_LIST =
|
||||
{
|
||||
{"is_one_works", is_one_works },
|
||||
{"roots_of_unity_is_the_expected_size", roots_of_unity_is_the_expected_size},
|
||||
{"roots_of_unity_are_plausible", roots_of_unity_are_plausible},
|
||||
{"reverse_works", reverse_works},
|
||||
{"expand_roots_is_plausible", expand_roots_is_plausible},
|
||||
{"new_fft_settings_is_plausible", new_fft_settings_is_plausible},
|
||||
{"is_power_of_two_works", is_power_of_two_works},
|
||||
{"fr_from_uint64_works", fr_from_uint64_works},
|
||||
{"fr_equal_works", fr_equal_works},
|
||||
{"compare_sft_fft", compare_sft_fft},
|
||||
{"roundtrip_fft", roundtrip_fft},
|
||||
{"inverse_fft", inverse_fft},
|
||||
|
73
src/fft_util.c
Normal file
73
src/fft_util.c
Normal file
@ -0,0 +1,73 @@
|
||||
/*
|
||||
* 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 "fft_util.h"
|
||||
|
||||
bool is_one(const blst_fr *fr_p) {
|
||||
uint64_t a[4];
|
||||
blst_uint64_from_fr(a, fr_p);
|
||||
return a[0] == 1 && a[1] == 0 && a[2] == 0 && a[3] == 0;
|
||||
}
|
||||
|
||||
bool is_power_of_two(uint64_t n) {
|
||||
return (n & (n - 1)) == 0;
|
||||
}
|
||||
|
||||
void fr_from_uint64(blst_fr *a, uint64_t n) {
|
||||
uint64_t vals[] = {n, 0, 0, 0};
|
||||
blst_fr_from_uint64(a, vals);
|
||||
}
|
||||
|
||||
// Returns an array of powers of the root of unity
|
||||
// Allocates space for the array that needs to be freed later
|
||||
blst_fr *expand_root_of_unity(blst_fr *root_of_unity, uint64_t width) {
|
||||
blst_fr *roots = malloc((width + 1) * sizeof(blst_fr));
|
||||
roots[0] = one;
|
||||
roots[1] = *root_of_unity;
|
||||
|
||||
for (int i = 2; !is_one(&roots[i - 1]); i++) {
|
||||
assert(i <= width);
|
||||
blst_fr_mul(&roots[i], &roots[i - 1], root_of_unity);
|
||||
}
|
||||
assert(is_one(&roots[width]));
|
||||
|
||||
return roots;
|
||||
}
|
||||
|
||||
// Return a reversed copy of the list of Fr provided
|
||||
// `width` is one less than the length of `r`
|
||||
// Allocates space for the array that needs to be freed later
|
||||
blst_fr *reverse(blst_fr *r, uint64_t width) {
|
||||
blst_fr *rr = malloc((width + 1) * sizeof(blst_fr));
|
||||
for (int i = 0; i <= width; i++) {
|
||||
rr[i] = r[width - i];
|
||||
}
|
||||
return rr;
|
||||
}
|
||||
|
||||
FFTSettings new_fft_settings(unsigned int max_scale) {
|
||||
FFTSettings s;
|
||||
s.max_width = (uint64_t)1 << max_scale;
|
||||
blst_fr_from_uint64(&s.root_of_unity, scale2_root_of_unity[max_scale]);
|
||||
s.expanded_roots_of_unity = expand_root_of_unity(&s.root_of_unity, s.max_width);
|
||||
s.reverse_roots_of_unity = reverse(s.expanded_roots_of_unity, s.max_width);
|
||||
return s;
|
||||
}
|
||||
|
||||
void free_fft_settings(FFTSettings *s) {
|
||||
free(s->expanded_roots_of_unity);
|
||||
free(s->reverse_roots_of_unity);
|
||||
}
|
78
src/fft_util.h
Normal file
78
src/fft_util.h
Normal file
@ -0,0 +1,78 @@
|
||||
/*
|
||||
* 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 <stdlib.h>
|
||||
#include "c-kzg.h"
|
||||
|
||||
// This is 1 in Blst's `blst_fr` limb representation. Crazy but true.
|
||||
static const blst_fr one =
|
||||
{0x00000001fffffffeL, 0x5884b7fa00034802L, 0x998c4fefecbc4ff5L, 0x1824b159acc5056fL};
|
||||
|
||||
// MODULUS = 52435875175126190479447740508185965837690552500527637822603658699938581184513
|
||||
// PRIMITIVE_ROOT = 5
|
||||
// [pow(PRIMITIVE_ROOT, (MODULUS - 1) // (2**i), MODULUS) for i in range(32)]
|
||||
//
|
||||
// These are not in `blst_fr` limb format and must be converted via `blst_fr_from_uint64()`
|
||||
static const uint64_t scale2_root_of_unity[][4] =
|
||||
{
|
||||
{0x0000000000000001L, 0x0000000000000000L, 0x0000000000000000L, 0x0000000000000000L},
|
||||
{0xffffffff00000000L, 0x53bda402fffe5bfeL, 0x3339d80809a1d805L, 0x73eda753299d7d48L},
|
||||
{0x0001000000000000L, 0xec03000276030000L, 0x8d51ccce760304d0L, 0x0000000000000000L},
|
||||
{0x8dd702cb688bc087L, 0xa032824078eaa4feL, 0xa733b23a98ca5b22L, 0x3f96405d25a31660L},
|
||||
{0x95df4b360411fe73L, 0x1ef4e672ebc1e1bbL, 0x6e92a9c460afca4aL, 0x4f2c596e753e4fccL},
|
||||
{0x9733f3a6ba60eaa6L, 0xbd7fdf9c77487ae7L, 0xd84f8612c8b6cc00L, 0x476fa2fb6162ffabL},
|
||||
{0xd2fc5e69ac5db47fL, 0xa12a70a615d0b8e4L, 0x293b1d67bc8de5d9L, 0x0e4840ac57f86f5eL},
|
||||
{0xb750da4cab28e208L, 0x501dff643be95635L, 0x8cbe2437f0b4b276L, 0x07d0c802a94a946eL},
|
||||
{0x2cabadec2fe322b8L, 0x752c84f315412560L, 0x32a732ae1a3b0aefL, 0x2e95da59a33dcbf2L},
|
||||
{0x33811ea1fe0c65f4L, 0x15c1ad4c687f28a2L, 0xecfbede342dee7f4L, 0x1bb466679a5d88b1L},
|
||||
{0xd58a5af42d010ff9L, 0x79efd6b0570bf109L, 0x3ed6d55a6350721dL, 0x2f27b09858f43cefL},
|
||||
{0x74a1f6718c130477L, 0xa534af14b61e0abeL, 0xeb674a1a620890d7L, 0x43527a8bca252472L},
|
||||
{0x450d9f977ea8ee05L, 0x565af17137d56fc0L, 0xe155cb4893f9e9acL, 0x110cebd0c8e9101bL},
|
||||
{0x23c9159959a0be92L, 0x87d188ce7a027759L, 0x70491431cab3c3ccL, 0x0ac00eb8b3f7f8daL},
|
||||
{0x13e96ade69583404L, 0x82c057275306243dL, 0x77e48bf529ca9f2aL, 0x50646ac81fe19595L},
|
||||
{0xe6a354dda97eccd4L, 0x39929d2e88fbbc57L, 0xa22ba63dd6e7b1c8L, 0x42c22911f5f07f43L},
|
||||
{0x137b458acfc35f7aL, 0x0caba63a29c01b06L, 0x0409ee987a02402cL, 0x6709c6cd56aa725bL},
|
||||
{0x10251f7d8831e03eL, 0x77d85a937ff858ecL, 0xebe905bd4fb9ac5cL, 0x05deb333f8727901L},
|
||||
{0xbf87b689b9009408L, 0x4f730e7ddd3ccc96L, 0xfd7f05ba4610300cL, 0x5ef5e8db0b8ac903L},
|
||||
{0x6499688417cd0c14L, 0xa672867368812f7fL, 0x2e1d9a1922cc3253L, 0x3a689e83aa0a1d80L},
|
||||
{0x20b53cbe41144deaL, 0x870c46fac2f0fcbdL, 0x556c35f6537d6971L, 0x3436287f5f686d91L},
|
||||
{0x007e082a436ba2e7L, 0x67c6630f9116e877L, 0x36f8f165fb4460f7L, 0x6eee34d57e7046e0L},
|
||||
{0xc5b670eea53a56d1L, 0x127d1f4253037d7bL, 0x57d4257ea722c2e2L, 0x03ae26a333cbd838L},
|
||||
{0x1e91484876504cf8L, 0x55bbbf1eb63edd02L, 0xbcdafec84e55aa02L, 0x5145c4cd2dc0beb0L},
|
||||
{0x5b90153a1ab70e2cL, 0x8deffa3175fb0ab8L, 0xc553ae2346900c95L, 0x1d31dcdc6bd3118cL},
|
||||
{0x801c894c59a2e8ebL, 0xbc535c5ce12fc974L, 0x95508d2747d39803L, 0x16d9d3cdac5d094fL},
|
||||
{0x810fa372cca1d8beL, 0xc67b8c2882e0bfa7L, 0xdbb4edf0e2d35bc2L, 0x712d15805087c995L},
|
||||
{0xeb162203fd88f133L, 0xac96c38ff010ea74L, 0x4307987fe64cfc70L, 0x350fe98d37b7a114L},
|
||||
{0xaba2f51842f2a254L, 0x4d7f3c3aa71efc0cL, 0x97ae418dd274a80aL, 0x2967385d5e3e7682L},
|
||||
{0x75c55c7b575a0b79L, 0x3ba4a15774a7ded1L, 0xc3974d73a04fccf3L, 0x705aba4f4a939684L},
|
||||
{0x8409a9ea14ebb608L, 0xfad0084e66bac611L, 0x04287254811c1dfbL, 0x086d072b23b30c29L},
|
||||
{0xb427c9b367e4756aL, 0xc7537fb902ebc38dL, 0x51de21becd6a205fL, 0x6064ab727923597dL}
|
||||
};
|
||||
|
||||
typedef struct {
|
||||
uint64_t max_width;
|
||||
blst_fr root_of_unity;
|
||||
blst_fr *expanded_roots_of_unity;
|
||||
blst_fr *reverse_roots_of_unity;
|
||||
} FFTSettings;
|
||||
|
||||
bool is_one(const blst_fr *fr_p);
|
||||
bool is_power_of_two(uint64_t n);
|
||||
void fr_from_uint64(blst_fr *a, uint64_t n);
|
||||
blst_fr *expand_root_of_unity(blst_fr *root_of_unity, uint64_t width);
|
||||
blst_fr *reverse(blst_fr *r, uint64_t width);
|
||||
FFTSettings new_fft_settings(unsigned int max_scale);
|
||||
void free_fft_settings(FFTSettings *s);
|
146
src/fft_util_test.c
Normal file
146
src/fft_util_test.c
Normal file
@ -0,0 +1,146 @@
|
||||
/*
|
||||
* 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_util.h"
|
||||
|
||||
#define NUM_ROOTS 32
|
||||
|
||||
void is_one_works(void) {
|
||||
TEST_CHECK(true == is_one(&one));
|
||||
}
|
||||
|
||||
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_are_plausible(void) {
|
||||
blst_fr r;
|
||||
for (unsigned int i = 0; i < NUM_ROOTS; i++) {
|
||||
blst_fr_from_uint64(&r, scale2_root_of_unity[i]);
|
||||
for (unsigned int j = 0; j < i; j++) {
|
||||
blst_fr_sqr(&r, &r);
|
||||
}
|
||||
TEST_CHECK(true == is_one(&r));
|
||||
TEST_MSG("Root %d failed", i);
|
||||
}
|
||||
}
|
||||
|
||||
void reverse_works(void) {
|
||||
int n = 24;
|
||||
blst_fr arr[n + 1];
|
||||
blst_fr *rev, diff;
|
||||
|
||||
// Initialise - increasing values
|
||||
arr[0] = one;
|
||||
for (int i = 1; i <= n; i++) {
|
||||
blst_fr_add(arr + i, arr + i - 1, &one);
|
||||
}
|
||||
|
||||
// Reverse
|
||||
rev = reverse(arr, n);
|
||||
|
||||
// Verify - decreasing values
|
||||
for (int i = 0; i < n; i++) {
|
||||
blst_fr_sub(&diff, rev + i, rev + i + 1);
|
||||
TEST_CHECK(true == is_one(&diff));
|
||||
}
|
||||
TEST_CHECK(true == is_one(rev + n));
|
||||
|
||||
free(rev);
|
||||
}
|
||||
|
||||
void expand_roots_is_plausible(void) {
|
||||
// Just test one (largeish) value of scale
|
||||
unsigned int scale = 20;
|
||||
unsigned int width = 1 << scale;
|
||||
blst_fr root, *expanded, prod;
|
||||
|
||||
// Initialise
|
||||
blst_fr_from_uint64(&root, scale2_root_of_unity[scale]);
|
||||
expanded = expand_root_of_unity(&root, width);
|
||||
|
||||
// Verify - each pair should multiply to one
|
||||
TEST_CHECK(true == is_one(expanded + 0));
|
||||
TEST_CHECK(true == is_one(expanded + width));
|
||||
for (unsigned int i = 1; i <= width / 2; i++) {
|
||||
blst_fr_mul(&prod, expanded + i, expanded + width - i);
|
||||
TEST_CHECK(true == is_one(&prod));
|
||||
}
|
||||
|
||||
free(expanded);
|
||||
}
|
||||
|
||||
void new_fft_settings_is_plausible(void) {
|
||||
// Just test one (largeish) value of scale
|
||||
unsigned int scale = 21;
|
||||
unsigned int width = 1 << scale;
|
||||
blst_fr prod;
|
||||
FFTSettings s = new_fft_settings(scale);
|
||||
|
||||
// Verify - each pair should multiply to one
|
||||
for (unsigned int i = 1; i <= width; i++) {
|
||||
blst_fr_mul(&prod, s.expanded_roots_of_unity + i, s.reverse_roots_of_unity + i);
|
||||
TEST_CHECK(true == is_one(&prod));
|
||||
}
|
||||
|
||||
free_fft_settings(&s);
|
||||
}
|
||||
|
||||
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 fr_from_uint64_works(void) {
|
||||
blst_fr a;
|
||||
fr_from_uint64(&a, 1);
|
||||
TEST_CHECK(true == is_one(&a));
|
||||
}
|
||||
|
||||
void fr_equal_works(void) {
|
||||
blst_fr a, b;
|
||||
blst_fr_from_uint64(&a, scale2_root_of_unity[15]);
|
||||
blst_fr_from_uint64(&b, scale2_root_of_unity[16]);
|
||||
TEST_CHECK(true == fr_equal(&a, &a));
|
||||
TEST_CHECK(false == fr_equal(&a, &b));
|
||||
}
|
||||
|
||||
TEST_LIST =
|
||||
{
|
||||
{"is_one_works", is_one_works },
|
||||
{"roots_of_unity_is_the_expected_size", roots_of_unity_is_the_expected_size},
|
||||
{"roots_of_unity_are_plausible", roots_of_unity_are_plausible},
|
||||
{"reverse_works", reverse_works},
|
||||
{"expand_roots_is_plausible", expand_roots_is_plausible},
|
||||
{"new_fft_settings_is_plausible", new_fft_settings_is_plausible},
|
||||
{"is_power_of_two_works", is_power_of_two_works},
|
||||
{"fr_from_uint64_works", fr_from_uint64_works},
|
||||
{"fr_equal_works", fr_equal_works},
|
||||
{ NULL, NULL } /* zero record marks the end of the list */
|
||||
};
|
Loading…
x
Reference in New Issue
Block a user