test/test_hist.c

   1 /* Apache 2.0 INS-AMU 2015 */
   2
   3 #include <stdlib.h>
   4 #include <string.h>
   5
   6 #include "sddekit.h"
   7
   8 static void hist_t_fill(void *data, int n, double *t, int *indices, double *buf) {
   9         /* suppress unused arguments */
  10         (void) data; (void) indices;
  11         /* fill buffer from times */
  12         memcpy(buf, t, n*sizeof(double));
  13 }
  14
  15 #define ND 4
  16
  17 TEST(hist, basic) {
  18         int i, vi[ND];
  19         double dt, vd[ND], x[ND];
  20         sk_hist *h;
  21
  22         /* set up */
  23         dt = 0.1;
  24         vi[0] = 1;
  25         vi[1] = 1;
  26         vi[2] = 0;
  27         vi[3] = 0;
  28         vd[0] = 5.5 * dt;
  29         vd[1] = 4.5 * dt;
  30         vd[2] = 33.3 * dt;
  31         vd[3] = 0.0;
  32         h = sk_hist_alloc();
  33
  34         sk_hist_init(h, ND, vi, vd, 0.0, dt);
  35         EXPECT_EQ(ND, sk_hist_get_nd(h));
  36         EXPECT_EQ(0.0, sk_hist_get_t(h));
  37         EXPECT_EQ(dt, sk_hist_get_dt(h));
  38         EXPECT_EQ(2, sk_hist_get_nu(h));
  39         EXPECT_EQ(0,sk_hist_get_lim(h, 0));
  40         EXPECT_EQ(36,sk_hist_get_lim(h, 1));
  41         EXPECT_EQ(8 + 36,sk_hist_get_lim(h, 2));
  42         EXPECT_EQ(36, sk_hist_get_len(h, 0));
  43         EXPECT_EQ(8, sk_hist_get_len(h, 1));
  44         EXPECT_EQ(0, sk_hist_get_pos(h, 0));
  45         EXPECT_EQ(0, sk_hist_get_pos(h, 1));
  46         EXPECT_EQ(0, sk_hist_get_uvi(h, 0));
  47         EXPECT_EQ(1, sk_hist_get_uvi(h, 1));
  48         EXPECT_EQ(1, sk_hist_get_maxvi(h));
  49         EXPECT_EQ(0, sk_hist_get_vi2i(h, 0));
  50         EXPECT_EQ(1, sk_hist_get_vi2i(h, 1));
  51         EXPECT_EQ(vd[0], sk_hist_get_maxd(h, 1));
  52         EXPECT_EQ(vd[2], sk_hist_get_maxd(h, 0));
  53         EXPECT_EQ(1, sk_hist_get_vi(h, 0));
  54         EXPECT_EQ(1, sk_hist_get_vi(h, 1));
  55         EXPECT_EQ(0, sk_hist_get_vi(h, 2));
  56         EXPECT_EQ(0, sk_hist_get_vi(h, 3));
  57         EXPECT_EQ(vd[0], sk_hist_get_vd(h, 0));
  58         EXPECT_EQ(vd[1], sk_hist_get_vd(h, 1));
  59         EXPECT_EQ(vd[2], sk_hist_get_vd(h, 2));
  60         EXPECT_EQ(vd[3], sk_hist_get_vd(h, 3));
  61         EXPECT_TRUE(!sk_hist_buf_is_null(h));
  62
  63         sk_hist_fill(h, &hist_t_fill, NULL);
  64         for (i=0; i<35; i++)
  65                 EXPECT_EQ(-i*dt, sk_hist_get_buf_lin(h, i));
  66         EXPECT_EQ(dt, sk_hist_get_buf_lin(h, 35));
  67         for (i=0; i<7; i++)
  68                 EXPECT_EQ(-i*dt, sk_hist_get_buf_lin(h, 36+i));
  69         EXPECT_EQ(dt, sk_hist_get_buf_lin(h, 36+7));
  70
  71         sk_hist_get(h, dt/3, x);
  72         ASSERT_NEAR(x[0], -vd[0]+dt/3, 1e-15);
  73         ASSERT_NEAR(x[1], -vd[1]+dt/3, 1e-15);
  74         ASSERT_NEAR(x[2], -vd[2]+dt/3, 1e-15);
  75         ASSERT_NEAR(x[3], -vd[3]+dt/3, 1e-15);
  76
  77         x[0] = 1.5;
  78         x[1] = 1.0;
  79         sk_hist_set(h, dt/2.0, x);
  80         ASSERT_NEAR( sk_hist_get_buf_lin(h, 35), 3.0, 1e-15);
  81         ASSERT_NEAR( sk_hist_get_buf_lin(h, 36 + 7), 2.0, 1e-15);
  82
  83         sk_hist_free(h);
  84 }
  85
  86 /* port of TVB's history test */
  87 static void hist_fill_ones(void *data, int n, double *t, int *indices, double *buf) {
  88         int i;
  89         /* suppress unused arguments */
  90         (void) data; (void) indices; (void) t;
  91         /* fill buffer from times */
  92         for (i=0; i<n; i++)
  93                 buf[i] = 1.0;
  94 }
  95
  96 static SK_DEFSYS(hist_exact_sys) {
  97         (void) data; (void) nx; (void) nc; (void) t; (void) Jf; (void) Jg;
  98         (void) Jce; (void) hist; (void) i;
  99         f[0] = x[0] + c[0];
 100         g[0] = 0.0;
 101         c[0] = x[0];
 102         return 0;
 103 }
 104
 105 /* The values are obtained by stepping through time according to the
 106  * connecting
 107  *
 108  * w = [ 0 1 0 0      d = [ .  1  .  .
 109  *       0 0 1 0            .  .  6  .
 110  *       0 0 0 1            .  .  . 11
 111  *       0 0 0 0 ]          .  .  .  . ]
 112  *
 113  *
 114  *    ...
 115  *    1  1  1  1
 116  * [  1  1  1  1    t = 0
 117  *    2  2  2  1        1
 118  *    3  3  3  1        2
 119  *    5  4  4  1        3
 120  *    8  5  5  1        4
 121  *   12  6  6  1        5
 122  *   17  7  7  1        6
 123  *   23  8  8  1        7
 124  *   30 10  9  1        8
 125  *   38 13 10  1        9
 126  *   ....
 127  *
 128  */
 129
 130 static double xexpect[12][4] = {
 131         {1, 1, 1, 1},
 132         {2, 2, 2, 1},
 133         {3, 3, 3, 1},
 134         {5, 4, 4, 1},
 135         {8, 5, 5, 1},
 136         {12, 6, 6, 1},
 137         {17, 7, 7, 1},
 138         {23, 8, 8, 1},
 139         {30, 10, 9, 1},
 140         {38, 13, 10, 1},
 141         {48, 17, 11, 1},
 142         {61, 22, 12, 1}
 143 };
 144
 145 static SK_DEFOUT(hist_exact_out) {
 146         int i;
 147         /* unused */ (void) nc; (void) c; (void) data;
 148         EXPECT_EQ(4,nx);
 149         for (i=0; i<4; i++)
 150                 EXPECT_EQ(x[i], xexpect[(int) t][i]);
 151         return t < 11;
 152 }
 153
 154 TEST(hist, exact) {
 155         int i, n=4, nnz, *Or, *Ic;
 156         double w[16], d[16], *nzw, *nzd, x0[4];
 157         sk_net_data *net;
 158         sk_solv *sol;
 159         sk_sch_id_data *schd;
 160
 161         schd = sk_sch_id_alloc();
 162         sol = sk_solv_alloc();
 163         x0[0] = 1;
 164         x0[1] = 1;
 165         x0[2] = 1;
 166         x0[3] = 1;
 167         for (i=0; i<16; i++) {
 168                 w[i] = 0.0;
 169                 d[i] = (double) i;
 170         }
 171         for (i=0; i<n-1; i++)
 172                 w[i*n+i+1] = 1;
 173
 174         /* 3 nnz, so hist expects c size 3, but net setups for c size 4 */
 175
 176         /* setup network */
 177         sk_sparse_from_dense(n, n, w, d, 0.0, &nnz, &Or, &Ic, &nzw, &nzd);
 178         EXPECT_EQ(0,Or[0]);
 179         EXPECT_EQ(1,Or[1]);
 180         EXPECT_EQ(2,Or[2]);
 181         EXPECT_EQ(3,Or[3]);
 182         EXPECT_EQ(1,Ic[0]);
 183         EXPECT_EQ(2,Ic[1]);
 184         EXPECT_EQ(3,Ic[2]);
 185         EXPECT_EQ(1,nzw[0]);
 186         EXPECT_EQ(1,nzw[1]);
 187         EXPECT_EQ(1,nzw[2]);
 188         EXPECT_EQ(1,nzd[0]);
 189         EXPECT_EQ(6,nzd[1]);
 190         EXPECT_EQ(11,nzd[2]);
 191
 192         net = sk_net_alloc();
 193         sk_net_init1(net, n, hist_exact_sys, NULL, 1, 1, nnz, Or, Ic, nzw, nzd);
 194
 195         EXPECT_EQ(3,nnz);
 196         EXPECT_EQ(4,sk_net_get_ne(net));
 197
 198         /* setup scheme & driver */
 199         sk_sch_id_init(schd, n);
 200         sk_solv_init(sol, sk_net_sys, net, sk_sch_id, schd,
 201                         hist_exact_out, NULL, hist_fill_ones, NULL, 42,
 202                         n, x0, nnz, Ic, nzd, 0.0, 1.0);
 203
 204         EXPECT_EQ(3, sk_solv_get_nc(sol));
 205         EXPECT_EQ(3, sk_hist_get_nu(sk_solv_get_hist(sol)));
 206
 207         /* run */
 208         sk_solv_cont(sol);
 209
 210         /* tests performed in out function.. */
 211
 212         /* clean up */
 213         sk_solv_free(sol);
 214         sk_sch_id_free(schd);
 215         sk_net_free(net);
 216         sk_free(Or);
 217         sk_free(Ic);
 218         sk_free(nzw);
 219         sk_free(nzd);
 220 }