1 /*
2 Copyright (C) 2014 Abhinav Baid
3
4 This file is part of FLINT.
5
6 FLINT is free software: you can redistribute it and/or modify it under
7 the terms of the GNU Lesser General Public License (LGPL) as published
8 by the Free Software Foundation; either version 2.1 of the License, or
9 (at your option) any later version. See <https://www.gnu.org/licenses/>.
10 */
11
12 #include "fmpz_mat.h"
13 #include "ulong_extras.h"
14
15 #define FMPZ_MAT_CHOL_EPS (1.0E-9)
16
17 int
main(void)18 main(void)
19 {
20 int i;
21 FLINT_TEST_INIT(state);
22
23 flint_printf("chol_d....");
24 fflush(stdout);
25
26 /* check RR^T = A */
27 for (i = 0; i < 1000 * flint_test_multiplier(); i++)
28 {
29 fmpz_mat_t A;
30 d_mat_t R, Rt, Atmp, Btmp;
31
32 slong m;
33
34 m = n_randint(state, 10);
35
36 fmpz_mat_init(A, m, m);
37
38 d_mat_init(R, m, m);
39 d_mat_init(Rt, m, m);
40 d_mat_init(Atmp, m, m);
41 d_mat_init(Btmp, m, m);
42
43 fmpz_mat_randtest(A, state, 10);
44 fmpz_mat_gram(A, A);
45 fmpz_mat_get_d_mat(Atmp, A);
46 d_mat_zero(R);
47
48 fmpz_mat_chol_d(R, A);
49 d_mat_transpose(Rt, R);
50
51 d_mat_mul_classical(Btmp, R, Rt);
52
53 if (!d_mat_approx_equal(Atmp, Btmp, FMPZ_MAT_CHOL_EPS))
54 {
55 flint_printf("FAIL:\n");
56 flint_printf("A:\n");
57 fmpz_mat_print_pretty(A);
58 flint_printf("R:\n");
59 d_mat_print(R);
60 flint_printf("R^T:\n");
61 d_mat_print(Rt);
62 flint_printf("Btmp:\n");
63 d_mat_print(Btmp);
64 abort();
65 }
66
67 fmpz_mat_clear(A);
68
69 d_mat_clear(R);
70 d_mat_clear(Rt);
71 d_mat_clear(Atmp);
72 d_mat_clear(Btmp);
73 }
74
75 FLINT_TEST_CLEANUP(state);
76
77 flint_printf("PASS\n");
78 return EXIT_SUCCESS;
79 }
80