half.c

   1 /*
   2  * actually do mapping from metaspace to realspace
   3  */
   4
   5 #include "half.h"
   6
   7 /*
   8  * Our metaspace dimensions are, in order:
   9  *   - y realspace coordinate of symmetry axis vertices
  10  *   - x,y,z realspace coordinates of even-y opposite-of-axis vertices
  11  *   - x,y,z realspace coordinates of even-y non-axis vertices
  12  *   - x,y,z realspace coordinates of odd-y vertices in canonical order
  13  */
  14
  15 #define FOR_MAP_DIMENSION__VERTEX_PAIR(body,mirror)     \
  16   k=0; sign=-1; { body } { mirror } i++;                \
  17   k=1; sign=+1; { body } { mirror } i++;                \
  18   k=2; sign=-1; { body } { mirror } i++;
  19
  20 /*#define FOR_MAP_DIMENSION_VPRIME              \
  21   int vprime= N-1-v
  22   */
  23
  24 #define FOR_MAP_DIMENSION(body, mirror, degenerate)             \
  25   do{                                                           \
  26     int i,v,vprime,k,sign; /* public */                         \
  27     int fmd_i, fmd_x, fmd_y;                                    \
  28     i= 0;                                                       \
  29     for (fmd_i=0; fmd_i<HALF_DEGRADED_NVERTICES; fmd_i++) {     \
  30       v= fmd_i*X*2;                                             \
  31       k=0;          { degenerate }                              \
  32       k=1; sign=+1; { body } i++;                               \
  33       k=2;          { degenerate }                              \
  34     }                                                           \
  35     for (fmd_y=0; fmd_y<Y; fmd_y+=2) {                          \
  36       v= fmd_y << YSHIFT | (X/2);                               \
  37       FOR_MAP_DIMENSION__VERTEX_PAIR(body,)                     \
  38     }                                                           \
  39     for (fmd_y=0; fmd_y<Y; fmd_y+=2)                            \
  40       for (fmd_x=1; fmd_x<X/2; fmd_x++) {                       \
  41         v= fmd_y << YSHIFT | fmd_x;                             \
  42         vprime= ((Y-1) << YSHIFT | X) - v;                      \
  43         FOR_MAP_DIMENSION__VERTEX_PAIR(body,mirror)             \
  44       }                                                         \
  45     for (fmd_y=1; fmd_y<Y; fmd_y+=2)                            \
  46       for (fmd_x=0; fmd_x<X/2; fmd_x++) {                       \
  47         v= fmd_y << YSHIFT | fmd_x;                             \
  48         vprime= ((Y-1) << YSHIFT | (X-1)) - v;                  \
  49         FOR_MAP_DIMENSION__VERTEX_PAIR(body,mirror)             \
  50       }                                                         \
  51     assert(i == DIM);                                           \
  52   }while(0)
  53
  54 void pmap_dimensions(const struct Vertices *vs) {
  55   int ixfrom[N][D3], count[N][D3], bad;
  56   int vc,kc;
  57
  58   bad= 0;
  59
  60   FOR_VERTEX(vc,INNER)
  61     FOR_COORD(kc) {
  62       count[vc][kc]= 0;
  63       ixfrom[vc][kc]= -1;
  64     }
  65
  66 #define CHECK_V(v) assert(v >= 0); assert(v < N);
  67 #define CHECK_NEAR(a,b)                                 \
  68   if (abs(a-b) >= 1e-6) { printf("    BAD!"); bad++; }  \
  69   putchar('\n');
  70
  71   FOR_MAP_DIMENSION(
  72                     (printf("i=%2d  v=%03x k=%d  %# 20g\n",
  73                             i,v,k,
  74                             vs->a[v][k]
  75                             ));
  76                     CHECK_V(v)
  77                     ixfrom[v][k]=i;
  78                     count[v][k]++;,
  79
  80 //                  FOR_MAP_DIMENSION_VPRIME;
  81                     (printf("i=%2d v'=%03x k=%d  %# 20g := %# 20g",
  82                             i,vprime,k,
  83                             vs->a[vprime][k], sign*vs->a[v][k]
  84                             ));
  85                     CHECK_V(vprime)
  86                     CHECK_NEAR(vs->a[vprime][k], sign*vs->a[v][k])
  87                     ixfrom[vprime][k]=-i;
  88                     count[vprime][k]++;,
  89
  90                     (printf("zero  v=%03x k=%d  %# 20g := 0",
  91                             v,k,
  92                             vs->a[v][k]
  93                             ));
  94                     CHECK_V(v)
  95                     CHECK_NEAR(vs->a[v][k], 0)
  96                     count[v][k]++;
  97                     );
  98
  99   FOR_VERTEX(vc,INNER)
 100     FOR_COORD(kc) {
 101       if (count[vc][kc]==1) continue;
 102       printf("BAD %03x %d count=%d ixfrom=%d\n",
 103              vc,kc, count[vc][kc], ixfrom[vc][kc]);
 104       bad++;
 105     }
 106   assert(!bad);
 107 }
 108
 109 void map_dimensions(const double metavertex[DIM], Vertices out) {
 110   FOR_MAP_DIMENSION(out[v     ][k]=      metavertex[i];,
 111
 112 //                  FOR_MAP_DIMENSION_VPRIME;
 113                     out[vprime][k]= sign*metavertex[i];,
 114
 115                     out[v     ][k]=      0;
 116                     );
 117 }
 118
 119 void unmap_dimensions(double metavertex[DIM], const struct Vertices *in) {
 120   FOR_MAP_DIMENSION(
 121                     metavertex[i]= in->a[v][k];,
 122                     ,
 123                     );
 124 }
 125
 126
 127 /*
 128  *  OLD NUMBERING
 129  *
 130  * Our numbering is as follows:
 131  *
 132  *          :axis of symmetry
 133  *        | :
 134  *        | :
 135  *     ___| 0  ___  1  ___  2  ___  3  ___  4  ___    up to X"-1
 136  *        | :                                           where X"=(X/2)
 137  *        / :\    /  \    /  \    /  \    /  \
 138  *       /| : \  /    \  /    \  /    \  /    \
 139  *        |___ 0   ___ 1   ___ 2   ___ 3   ___ 4
 140  *        | :  +X"     +X"     +X"     +X"     +X"
 141  *       \| : /  \    /  \    /  \    /  \    /
 142  *        \ :/    \  /    \  /    \  /    \  /
 143  *     ___|2*X"___ 1+  ___ 2+  ___ 3+  __  4+  ___
 144  *        | :     2*X"    2*X"    2*X"    2*X"
 145  *        / :\    /  \    /  \    /  \    /  \
 146  *       /| : \  /    \  /    \  /    \  /    \
 147  *        |___ 0+  ___ 1+  ___ 2+  ___ 3+  ___ 4+
 148  *        | : 3*X"    3*X"    3*X"    3*X"    3*X"
 149  *        | :
 150  *
 151  * X" is still a power of two so a vertex pair can be named by a
 152  * number structured like this:
 153  *
 154  * vertex number":   0000 | y     | x
 155  *                         YBITS   XBITS-1
 156  *
 157  */