Commit 4825fd5d authored by Matthieu Schaller's avatar Matthieu Schaller

Added unit test to verify that the random numbers are generated uniformly over...

Added unit test to verify that the random numbers are generated uniformly over [0,1] for all IDs and time-steps.
parent 3ebdebec
......@@ -110,6 +110,7 @@ tests/testInteractions
tests/testInteractions.sh
tests/testSymmetry
tests/testMaths
tests/testRandom
tests/testThreadpool
tests/testParser
tests/parser_output.yml
......
......@@ -66,7 +66,7 @@ INLINE static double random_unit_interval(const long long int id,
with either the current integer time or the particle IDs.
The calculation overflows on purpose. */
unsigned int seed = ((937LL * id + 1109LL) % 2147987LL +
(ti_current - 1) % 1514917LL + (long long)type) %
(ti_current - 1LL) % 1514917LL + (long long)type) %
seed_range;
/* Generate a random number between 0 and 1. */
......
......@@ -61,6 +61,7 @@
#include "potential.h"
#include "profiler.h"
#include "queue.h"
#include "random.h"
#include "restart.h"
#include "runner.h"
#include "scheduler.h"
......
......@@ -23,7 +23,7 @@ AM_LDFLAGS = ../src/.libs/libswiftsim.a $(HDF5_LDFLAGS) $(HDF5_LIBS) $(FFTW_LIBS
TESTS = testGreetings testMaths testReading.sh testSingle testKernel \
testActivePair.sh test27cells.sh test27cellsPerturbed.sh \
testParser.sh testSPHStep test125cells.sh test125cellsPerturbed.sh testFFT \
testAdiabaticIndex \
testAdiabaticIndex testRandom \
testMatrixInversion testThreadpool testDump testLogger testInteractions.sh \
testVoronoi1D testVoronoi2D testVoronoi3D testGravityDerivatives \
testPeriodicBC.sh testPeriodicBCPerturbed.sh testPotentialSelf \
......@@ -34,7 +34,7 @@ TESTS = testGreetings testMaths testReading.sh testSingle testKernel \
# List of test programs to compile
check_PROGRAMS = testGreetings testReading testSingle testTimeIntegration \
testSPHStep testActivePair test27cells test27cells_subset test125cells testParser \
testKernel testFFT testInteractions testMaths \
testKernel testFFT testInteractions testMaths testRandom \
testSymmetry testThreadpool \
testAdiabaticIndex testRiemannExact testRiemannTRRS \
testRiemannHLLC testMatrixInversion testDump testLogger \
......@@ -51,6 +51,8 @@ testGreetings_SOURCES = testGreetings.c
testMaths_SOURCES = testMaths.c
testRandom_SOURCES = testRandom.c
testReading_SOURCES = testReading.c
testSelectOutput_SOURCES = testSelectOutput.c
......
/*******************************************************************************
* This file is part of SWIFT.
* Copyright (C) 2019 Matthieu Schaller (schaller@strw.leidenuniv.nl)
*
* This program is free software: you can redistribute it and/or modify
* it under the terms of the GNU Lesser General Public License as published
* by the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU Lesser General Public License
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*
******************************************************************************/
/* Config parameters. */
#include "../config.h"
#include <fenv.h>
/* Local headers. */
#include "swift.h"
int main(int argc, char* argv[]) {
/* Initialize CPU frequency, this also starts time. */
unsigned long long cpufreq = 0;
clocks_set_cpufreq(cpufreq);
/* Choke on FPEs */
#ifdef HAVE_FE_ENABLE_EXCEPT
feenableexcept(FE_DIVBYZERO | FE_INVALID | FE_OVERFLOW);
#endif
/* Get some randomness going */
const int seed = time(NULL);
message("Seed = %d", seed);
srand(seed);
/* Time-step size */
const int time_bin = 29;
/* Try a few different values for the ID */
for (int i = 0; i < 20; ++i) {
const long long id = rand() * (1LL << 31) + rand();
const integertime_t increment = (1LL << time_bin);
message("Testing id=%lld time_bin=%d", id, time_bin);
double total = 0., total2 = 0.;
int count = 0;
/* Check that the numbers are uniform over the full-range of useful
* time-steps */
for (integertime_t ti_current = 0LL; ti_current < max_nr_timesteps;
ti_current += increment) {
ti_current += increment;
const double r =
random_unit_interval(id, ti_current, random_number_star_formation);
total += r;
total2 += r * r;
count++;
}
const double mean = total / (double)count;
const double var = total2 / (double)count - mean * mean;
/* Verify that the mean and variance match the expected values for a uniform
* distribution */
if ((fabs(mean - 0.5) / 0.5 > 1e-4) ||
(fabs(var - 1. / 12.) / (1. / 12.) > 1e-4)) {
message("Test failed!");
message("Result: count=%d mean=%f var=%f", count, mean, var);
message("Expected: count=%d mean=%f var=%f", count, 0.5f, 1. / 12.);
return 1;
}
}
return 0;
}
Markdown is supported
0%
or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment