Second review after Fredrik's comments

This commit is contained in:
Albin Ahlbäck 2021-05-18 01:41:46 +02:00
parent 7eab089f7c
commit ef6b9f2449
9 changed files with 140 additions and 45 deletions

4
arb.h
View file

@ -333,9 +333,7 @@ void arb_randtest(arb_t x, flint_rand_t state, slong prec, slong mag_bits);
void arb_randtest_special(arb_t x, flint_rand_t state, slong prec, slong mag_bits);
void arb_randtest_uniform(arb_t x, flint_rand_t state, slong prec);
void arb_randtest_uniform_exact(arb_t x, flint_rand_t state, slong prec);
void arb_urandom(arb_t x, flint_rand_t state, slong prec, arf_rnd_t rnd);
void arb_add_error_arf(arb_t x, const arf_t err);

View file

@ -88,17 +88,3 @@ arb_randtest_special(arb_t x, flint_rand_t state, slong prec, slong mag_bits)
}
}
void
arb_randtest_uniform(arb_t x, flint_rand_t state, slong prec)
{
arf_randtest_uniform(arb_midref(x), state, prec);
mag_randtest_uniform(arb_radref(x), state);
}
void
arb_randtest_uniform_exact(arb_t x, flint_rand_t state, slong prec)
{
arf_randtest_uniform(arb_midref(x), state, prec);
mag_zero(arb_radref(x));
}

84
arb/test/t-urandom.c Normal file
View file

@ -0,0 +1,84 @@
/*
Copyright (C) 2021 Fredrik Johansson
This file is part of Arb.
Arb is free software: you can redistribute it and/or modify it under
the terms of the GNU Lesser General Public License (LGPL) as published
by the Free Software Foundation; either version 2.1 of the License, or
(at your option) any later version. See <http://www.gnu.org/licenses/>.
*/
#include "arb.h"
#define N 10000
int main()
{
slong iter;
slong prec;
flint_rand_t state;
arf_rnd_t round;
arb_t rand[N];
arb_t m; /* mean */
arb_t s; /* variance */
arb_t mp;
arb_t sp;
flint_printf("urandom....");
fflush(stdout);
flint_randinit(state);
arb_init(m);
arb_init(s);
arb_init(mp);
arb_init(sp);
prec = 299;
round = ARF_RND_DOWN;
for (iter = 0; iter < N; iter++)
{
arb_init(rand[iter]);
arb_urandom(rand[iter], state, prec, round);
arb_add(m, m, rand[iter], prec);
}
arb_div_si(m, m, N, prec);
for (iter = 0; iter < N; iter++)
{
arb_t tmp;
arb_init(tmp);
arb_sub(tmp, rand[iter], m, prec);
arb_sqr(tmp, tmp, prec);
arb_add(s, s, tmp, prec);
arb_clear(tmp);
}
arb_div_si(s, s, N, prec);
/* one percent deviation */
arb_set_str(mp, "0.5 +/- 0.005", prec);
arb_set_str(sp, "0.083333 +/- 0.00083", prec);
if (!arb_contains(mp, m))
{
flint_printf("FAIL: mean\n\n");
flint_printf("m = "); arb_printd(m, 15); flint_printf("\n\n");
flint_abort();
}
if (!arb_contains(sp, s))
{
flint_printf("FAIL: variance\n\n");
flint_printf("s = "); arb_printd(s, 15); flint_printf("\n\n");
flint_abort();
}
for (iter = 0; iter < N; iter++) arb_clear(rand[iter]);
arb_clear(m);
arb_clear(s);
flint_randclear(state);
flint_cleanup();
flint_printf("PASS\n");
return EXIT_SUCCESS;
}

20
arb/urandom.c Normal file
View file

@ -0,0 +1,20 @@
/*
Copyright (C) 2021 Fredrik Johansson
This file is part of Arb.
Arb is free software: you can redistribute it and/or modify it under
the terms of the GNU Lesser General Public License (LGPL) as published
by the Free Software Foundation; either version 2.1 of the License, or
(at your option) any later version. See <http://www.gnu.org/licenses/>.
*/
#include "arb.h"
void
arb_urandom(arb_t x, flint_rand_t state, slong prec, arf_rnd_t rnd)
{
arf_urandom(arb_midref(x), state, prec, rnd);
mag_zero(arb_radref(x));
}

2
arf.h
View file

@ -797,7 +797,7 @@ void arf_randtest_not_zero(arf_t x, flint_rand_t state, slong bits, slong mag_bi
void arf_randtest_special(arf_t x, flint_rand_t state, slong bits, slong mag_bits);
void arf_randtest_uniform(arf_t x, flint_rand_t state, slong bits);
void arf_urandom(arf_t x, flint_rand_t state, slong bits, arf_rnd_t rnd);
#define MUL_MPFR_MIN_LIMBS 25
#define MUL_MPFR_MAX_LIMBS 10000

View file

@ -56,16 +56,3 @@ arf_randtest_special(arf_t x, flint_rand_t state, slong bits, slong mag_bits)
}
}
void
arf_randtest_uniform(arf_t x, flint_rand_t state, slong bits)
{
slong prec = bits;
prec += 128;
fmpz_t t;
fmpz_init(t);
fmpz_randtest_unsigned(t, state, prec);
arf_set_round_fmpz(x, t, bits, ARF_RND_DOWN);
fmpz_sub_si(ARF_EXPREF(x), ARF_EXPREF(x), prec);
fmpz_clear(t);
}

34
arf/urandom.c Normal file
View file

@ -0,0 +1,34 @@
/*
Copyright (C) 2021 Fredrik Johansson
This file is part of Arb.
Arb is free software: you can redistribute it and/or modify it under
the terms of the GNU Lesser General Public License (LGPL) as published
by the Free Software Foundation; either version 2.1 of the License, or
(at your option) any later version. See <http://www.gnu.org/licenses/>.
*/
#include "arf.h"
void
arf_urandom(arf_t x, flint_rand_t state, slong bits, arf_rnd_t rnd)
{
slong prec = bits;
fmpz_t n;
fmpz_t t;
prec += 128;
fmpz_init(n);
fmpz_one(n);
fmpz_mul_2exp(n, n, (ulong) prec);
fmpz_init(t);
fmpz_randm(t, state, n);
arf_set_round_fmpz(x, t, bits, rnd);
arf_mul_2exp_si(x, x, -prec);
fmpz_clear(t);
}

2
mag.h
View file

@ -508,8 +508,6 @@ void mag_randtest_special(mag_t x, flint_rand_t state, slong expbits);
void mag_randtest(mag_t x, flint_rand_t state, slong expbits);
void mag_randtest_uniform(mag_t x, flint_rand_t state);
void mag_fprint(FILE * file, const mag_t x);
void mag_fprintd(FILE * file, const mag_t x, slong d);

View file

@ -45,15 +45,3 @@ mag_randtest(mag_t x, flint_rand_t state, slong expbits)
mag_zero(x);
}
void
mag_randtest_uniform(mag_t x, flint_rand_t state)
{
slong prec = 64;
fmpz_t t;
fmpz_init(t);
fmpz_randtest_unsigned(t, state, prec);
mag_set_fmpz_lower(x, t);
fmpz_sub_si(MAG_EXPREF(x), MAG_EXPREF(x), prec);
fmpz_clear(t);
}