#include "profiler.h"
#include "nmod_mat.h"
#include "ulong_extras.h"
#include "thread_support.h"
#include "perm.h"
typedef struct
{
slong n;
ulong modulus;
int algorithm;
} mat_lu_t;
void sample(void * arg, ulong count)
{
mat_lu_t * params = (mat_lu_t *) arg;
int algorithm = params->algorithm;
slong * P;
nmod_mat_t A, LU;
ulong i;
flint_rand_t state;
flint_rand_init(state);
nmod_mat_init(A, params->n, params->n, params->modulus);
nmod_mat_init(LU, params->n, params->n, params->modulus);
nmod_mat_randfull(A, state);
P = _perm_init(params->n);
prof_start();
if (algorithm == 0)
for (i = 0; i < count; i++)
{
nmod_mat_set(LU, A);
nmod_mat_lu(P, LU, 0);
}
else if (algorithm == 1)
for (i = 0; i < count; i++)
{
nmod_mat_set(LU, A);
nmod_mat_lu_classical(P, LU, 0);
}
else if (algorithm == 2)
for (i = 0; i < count; i++)
{
nmod_mat_set(LU, A);
nmod_mat_lu_classical_delayed(P, LU, 0);
}
else
for (i = 0; i < count; i++)
{
nmod_mat_set(LU, A);
nmod_mat_lu_recursive(P, LU, 0);
}
prof_stop();
nmod_mat_clear(A);
nmod_mat_clear(LU);
_perm_clear(P);
flint_rand_clear(state);
}
slong bits_tab[] = { 5, 14, 15, 25, 30, 31, 32, 33, 60, 61, 62, 63, 64, 0 };
int main(void)
{
double max;
mat_lu_t params;
slong i, bits, n;
flint_printf("nmod_mat_lu:\n");
flint_set_num_threads(1);
flint_printf("threads = %wd\n", flint_get_num_threads());
for (i = 0; (bits = bits_tab[i]) != 0; i++)
{
for (n = 4; n <= 1000; n += n/4 + 1)
{
double min_default = 0, min_classical = 0, min_delayed = 0, min_recursive = 0;
params.n = n;
params.modulus = n_nextprime(UWORD(1) << (bits - 1), 0);
params.algorithm = 0;
prof_repeat(&min_default, &max, sample, ¶ms);
params.algorithm = 1;
prof_repeat(&min_classical, &max, sample, ¶ms);
params.algorithm = 2;
prof_repeat(&min_delayed, &max, sample, ¶ms);
params.algorithm = 3;
prof_repeat(&min_recursive, &max, sample, ¶ms);
flint_printf("b = %wd n = %wd default %.2f us classical %.2f us delayed %.2f us recursive %.2f us\n",
bits, n, min_default, min_classical, min_delayed, min_recursive);
}
}
return 0;
}