sysrap/tests/smath_test.cc

0001 // ~/o/sysrap/tests/smath_test.sh
0002
0003 #include <iostream>
0004 #include <iomanip>
0005
0006 #include "NPX.h"
0007 #include "ssys.h"
0008 #include "scuda.h"
0009 #include "squad.h"
0010 #include "smath.h"
0011
0012 struct smath_test
0013 {
0014     static int count_nibbles();
0015     static int rotateUz();
0016     static int erfcinvf();
0017     static int main();
0018 };
0019
0020
0021 int smath_test::count_nibbles()
0022 {
0023     typedef unsigned long long ULL ;
0024     static const int N = 21 ;
0025     ULL xx[N] ;
0026     int nn[N] ;
0027
0028     xx[ 0] = 0x0123456789abcdefull ; nn[ 0] = 15 ;
0029     xx[ 1] = 0x0023456789abcdefull ; nn[ 1] = 14 ;
0030     xx[ 2] = 0x0003456789abcdefull ; nn[ 2] = 13 ;
0031     xx[ 3] = 0x0000456789abcdefull ; nn[ 3] = 12 ;
0032     xx[ 4] = 0x0000056789abcdefull ; nn[ 4] = 11 ;
0033     xx[ 5] = 0x0000006789abcdefull ; nn[ 5] = 10 ;
0034     xx[ 6] = 0x0000000789abcdefull ; nn[ 6] =  9 ;
0035     xx[ 7] = 0x0000000089abcdefull ; nn[ 7] =  8 ;
0036     xx[ 8] = 0x0000000009abcdefull ; nn[ 8] =  7 ;
0037     xx[ 9] = 0x0000000000abcdefull ; nn[ 9] =  6 ;
0038     xx[10] = 0x00000000000bcdefull ; nn[10] =  5 ;
0039     xx[11] = 0x000000000000cdefull ; nn[11] =  4 ;
0040     xx[12] = 0x0000000000000defull ; nn[12] =  3 ;
0041     xx[13] = 0x00000000000000efull ; nn[13] =  2 ;
0042     xx[14] = 0x000000000000000full ; nn[14] =  1 ;
0043     xx[15] = 0x0000000000000000ull ; nn[15] =  0 ;
0044     xx[16] = 0x0000d00e000a000dull ; nn[16] =  4 ;
0045     xx[17] = 0x0000100000000000ull ; nn[17] =  1 ;
0046     xx[18] = 0xa123456789abcdefull ; nn[18] = 16 ;
0047     xx[19] = 0x1111111111111111ull ; nn[19] = 16 ;
0048     xx[20] = 0x0000000000000000ull ; nn[20] =  0 ;
0049
0050     for(int i=0 ; i < N ; i++)
0051     {
0052         ULL x = xx[i] ;
0053         int n = smath::count_nibbles(x) ;
0054         std::cout
0055             << " i " << std::setw(3)  << i
0056             << " x " << std::setw(16) << std::hex << x  << std::dec
0057             << " n " << std::setw(3)  << n
0058             << " nn[i] " << std::setw(3) << nn[i]
0059             << std::endl
0060             ;
0061     }
0062     return 0;
0063 }
0064
0065
0066 struct float3_x3
0067 {
0068     float3 u ;
0069     float3 d0 ;
0070     float3 d1 ;
0071 };
0072
0073
0074 /**
0075 rotateUz
0076 ---------
0077
0078 Consider mom is some direction, say +Z::
0079
0080    (0, 0, 1)
0081
0082 There is a circle of vectors that are perpendicular
0083 to that mom, all in the XY plane, and with dot product
0084 with that direction of zero::
0085
0086    ( cos(phi), sin(phi), 0 )    phi 0->2pi
0087
0088 **/
0089
0090 int smath_test::rotateUz()
0091 {
0092     float3 u = normalize(make_float3( 1.f, 0.f, -1.f ));
0093     std::cout << " u " << u << std::endl ;
0094
0095     static const int N = 16 ;
0096     std::vector<float3_x3> vv(N+1) ;
0097     for(int i=0 ; i <= N ; i++)
0098     {
0099         float phi = 2.f*M_PIf*float(i)/float(N) ;
0100         float3 d0 = make_float3( cos(phi), sin(phi), 0.f ) ;
0101         // d0: ring of vectors in XY plane, "around" the +Z direction
0102
0103         float3 d1(d0);
0104         smath::rotateUz(d1,u);
0105
0106         vv[i].u = u ;
0107         vv[i].d0 = d0 ;
0108         vv[i].d1 = d1 ;
0109
0110         // d1: rotated XY ring of vectors to point in direction u
0111         // So all the d1 are perpendicular to u
0112
0113         std::cout
0114             << std::setw(2) << i
0115             << " d0 " << d0
0116             << " d1 " << d1
0117             << " dot(d1,u)*1e6 " << dot(d1,u)*1e6
0118             << std::endl
0119             ;
0120     }
0121
0122     NP* a = NPX::ArrayFromVec<float,float3_x3>(vv,3,3);
0123     a->save("$FOLD/rotateUz.npy");
0124
0125     return 0;
0126 }
0127
0128 int smath_test::erfcinvf()
0129 {
0130     float SQRT2 = sqrtf(2.f) ;
0131
0132     int N = 100 ;
0133     for(int i=0 ; i < N ; i++)
0134     {
0135         float u2 = 2.f*float(i)/float(N-1) ;
0136         float v = -SQRT2*smath::erfcinvf(u2) ;
0137         std::cout
0138             << " i " << std::setw(5) << i
0139             << " u2 " << std::setw(10) << std::setprecision(5) << std::fixed << u2
0140             << " v  " << std::setw(10) << std::setprecision(5) << std::fixed << v
0141             << std::endl
0142             ;
0143     }
0144     return 0;
0145 }
0146
0147 int smath_test::main()
0148 {
0149     const char* TEST = ssys::getenvvar("TEST","rotateUz");
0150     bool ALL = strcmp(TEST, "ALL") == 0 ;
0151
0152     int rc = 0;
0153     if(ALL||0==strcmp(TEST,"count_nibbles")) rc += count_nibbles();
0154     if(ALL||0==strcmp(TEST,"rotateUz")) rc += rotateUz();
0155     if(ALL||0==strcmp(TEST,"erfcinvf")) rc += erfcinvf();
0156
0157     return rc;
0158 }
0159
0160 int main(int argc, char** argv)
0161 {
0162     return smath_test::main();
0163 }
0164 // ~/o/sysrap/tests/smath_test.sh
0165