test/linalg_cho.h

   1
   2 #define MAIN(x) linalg_cho##x
   3 #include "test.h"
   4 #include "a/linalg.h"
   5
   6 static void show(a_float const *p, a_size m, a_size n)
   7 {
   8     a_size i, ii;
   9     for (i = 0; i < m; ++i)
  10     {
  11         debug("\t");
  12         for (ii = 0; ii < n; ++ii)
  13         {
  14             a_float const x = *p++;
  15             if (ii) { debug(","); }
  16             debug("%+.6" A_FLOAT_PRI "f", x);
  17         }
  18         debug("\n");
  19     }
  20 }
  21
  22 int main(int argc, char *argv[]) /* NOLINT(misc-definitions-in-headers) */
  23 {
  24 #define N 3
  25     int i, ii;
  26     a_float b[N];
  27     a_float m[N * N];
  28     a_float a[N * N] = {4, 12, -16, 12, 37, -43, -16, -43, 98};
  29     main_init(argc, argv, 1);
  30     debug("A=\n");
  31     show(a, N, N);
  32     if (a_linalg_cho(a, N) == 0)
  33     {
  34         a_linalg_cho_get_L(a, N, m);
  35         debug("L=\n");
  36         show(m, N, N);
  37         a_linalg_cho_inv(a, N, b, m);
  38         debug("INV=\n");
  39         show(m, N, N);
  40         for (i = 0; i < N; ++i)
  41         {
  42             a_float *x = m + a_size_c(i) * N;
  43             for (ii = 0; ii < N; ++ii)
  44             {
  45                 x[ii] = 0;
  46             }
  47             x[i] = 1;
  48             a_linalg_cho_solve(a, N, x);
  49         }
  50         a_linalg_T1(m, N);
  51         debug("I=\n");
  52         show(m, N, N);
  53         a_linalg_T2(a, m, N, N);
  54         debug("T=\n");
  55         show(a, N, N);
  56         if (a_linalg_cho(m, N) == 0)
  57         {
  58             a_linalg_cho_get_L(m, N, a);
  59             debug("L=\n");
  60             show(a, N, N);
  61             a_linalg_cho_inv(m, N, b, a);
  62             debug("INV=\n");
  63             show(a, N, N);
  64         }
  65     }
  66     return 0;
  67 }