diff --git a/devmods/reps/counter.cpp b/devmods/reps/counter.cpp new file mode 100644 index 00000000..8480ba19 --- /dev/null +++ b/devmods/reps/counter.cpp @@ -0,0 +1,73 @@ +// a float-like type to count operations + +namespace reps { + +std::map counts; + +#define C(x) {} + +std::array cbc; + +constexpr int cbcAdd = 1; +constexpr int cbcMul = 2; +constexpr int cbcDiv = 3; +constexpr int cbcTrig = 4; + +struct countfloat { + ld x; + + countfloat() {} + explicit countfloat(ld _x) : x(_x) {} + explicit operator ld() { return x; } + operator bool() { return x != 0; } + + countfloat operator +(countfloat a) const { C("plus"); cbc[1]++; return countfloat(x + a.x); } + countfloat operator -(countfloat a) const { C("plus"); cbc[1]++; return countfloat(x + a.x); } + countfloat operator -() const { return countfloat(-x); } + countfloat operator +() const { return countfloat(+x); } + countfloat operator *(countfloat a) const { C("mul"); cbc[2]++; return countfloat(x * a.x); } + countfloat operator /(countfloat a) const { C("div"); cbc[3]++; return countfloat(x / a.x); } + + bool operator <(countfloat a) const { return x < a.x; } + bool operator >(countfloat a) const { return x > a.x; } + bool operator <=(countfloat a) const { return x <= a.x; } + bool operator >=(countfloat a) const { return x >= a.x; } + + countfloat& operator +=(countfloat a) { C("plus"); cbc[1]++; x += a.x; return self; } + countfloat& operator -=(countfloat a) { C("plus"); cbc[1]++; x -= a.x; return self; } + countfloat& operator *=(countfloat a) { C("mul"); cbc[2]++; x *= a.x; return self; } + countfloat& operator /=(countfloat a) { C("mul"); cbc[2]++; x /= a.x; return self; } + countfloat& operator *=(int a) { if(a != 1 && a != -1) C("mul"+hr::its(a)); x *= a; return self; } + countfloat& operator /=(int a) { if(a != 1 && a != -1) C("div"+hr::its(a)); x /= a; return self; } + + friend countfloat sin(countfloat a) { cbc[4]++; C("sin"); return countfloat(sin(a.x)); } + friend countfloat cos(countfloat a) { cbc[4]++; C("cos"); return countfloat(cos(a.x)); } + friend countfloat tan(countfloat a) { cbc[4]++; C("tan"); return countfloat(tan(a.x)); } + friend countfloat sinh(countfloat a) { cbc[4]++; C("sinh"); return countfloat(sinh(a.x)); } + friend countfloat cosh(countfloat a) { cbc[4]++; C("cosh"); return countfloat(cosh(a.x)); } + friend countfloat tanh(countfloat a) { cbc[4]++; C("cosh"); return countfloat(tanh(a.x)); } + friend countfloat asinh(countfloat a) { cbc[4]++; C("asinh"); return countfloat(asinh(a.x)); } + friend countfloat acosh(countfloat a) { cbc[4]++; C("acosh"); return countfloat(acosh(a.x)); } + friend countfloat acos(countfloat a) { cbc[4]++; C("acos"); return countfloat(acos(a.x)); } + friend countfloat exp(countfloat a) { cbc[4]++; C("exp"); return countfloat(exp(a.x)); } + friend countfloat log(countfloat a) { cbc[4]++; C("log"); return countfloat(log(a.x)); } + friend countfloat sqrt(countfloat a) { cbc[4]++; C("sqrt"); return countfloat(sqrt(a.x)); } + friend countfloat atan2(countfloat a, countfloat b) { cbc[4]++; C("atan"); return countfloat(atan2(a.x, b.x)); } + friend countfloat pow(countfloat a, ld b) { cbc[4]++; C("pow" + hr::fts(b)); return countfloat(pow(a.x, b)); } + friend countfloat abs(countfloat a) { return countfloat(abs(a.x)); } + countfloat operator *(int a) const { if(a != 1 && a != -1) C("mul" + hr::its(a)); return countfloat(x * a); } + countfloat operator /(int a) const { if(a != 1 && a != -1) C("div" + hr::its(a)); return countfloat(x / a); } + + friend bool isinf(countfloat a) { return isinf(a.x); } + friend bool isnan(countfloat a) { return isnan(a.x); } + }; + +template<> countfloat get_deg (int deg) { return countfloat( M_PI * deg / 180 ); } + +} + +namespace hr { + void print(hr::hstream& hs, ::reps::countfloat b) { + print(hs, b.x); + } +} diff --git a/devmods/reps/multivector.cpp b/devmods/reps/multivector.cpp new file mode 100644 index 00000000..20062cff --- /dev/null +++ b/devmods/reps/multivector.cpp @@ -0,0 +1,330 @@ +namespace reps { + +TD struct mvector { + array values; + typename D::Number& operator [] (int i) { return values[i]; } + const typename D::Number& operator [] (int i) const { return values[i]; } + mvector operator + (const mvector& M) const { + mvector result; + for(int i=0; i, D::Dim> values; + array& operator [] (int i) { return values[i]; } + const array& operator [] (int i) const { return values[i]; } + matrix operator * (const matrix& M) const { + matrix result; + for(int i=0; i operator * (const mvector& V) const { + mvector result; + for(int i=0; i zero_vector() { + mvector result; + for(int i=0; i unit_vector(int id) { + mvector result; + for(int i=0; i>; + +TD std::string nz(const multivector& a) { + constexpr int mdim = 1< Number(1e-10)) { + if(str.s != "") print(str, " "); + if(a[i] > Number(0)) print(str, "+"); + print(str, a[i]); + for(int u=0; u unit(const typename D::Number& a) { + auto res = zero_vector>(); + res[0] = a; + return res; + } + +TD constexpr multivector embed(const mvector& a) { + auto res = zero_vector>(); + for(int i=0; i unembed(const multivector& a) { + mvector res; + for(int i=0; i multimul(const multivector& a, const multivector& b) { + constexpr int mdim = 1<>(); + for(mvindex i=0; i(i, j); + } + return res; + } + +template +multivector chkmul(const multivector& a, const multivector& b) { + constexpr int mdim = 1<>(); + + /* we initialize with 0s and then add stuff, so one add per component is not necessary */ + for(mvindex i=0; i(i, j); + else if(B::isflat(j)) + res[i^j] += a[i] * mul_sign(i, j); + else + res[i^j] += a[i] * b[j] * mul_sign(i, j); + } + return res; + } + +TD multivector conjugate(const multivector& a) { + constexpr int mdim = 1<(i); + return res; + } + +TD multivector transpose(const multivector& a) { + constexpr int mdim = 1<(i); + return res; + } + +template multivector apply_nm(multivector a); + +TD using poincare_rotation = std::pair, multivector>; + +/** decompose o into the poincare part and the rotational component */ +TD poincare_rotation despin2(const multivector& a) { + auto p = a; + for(int i=(1<<(D::Dim-1)); i<(1<<(D::Dim)); i++) p[i] = typename D::Number(0); + p = p * pow(chkmul,rotational,units>(p, conjugate(p))[0], -0.5); + auto p1 = chkmul,rotational,poincare>(a, conjugate(p)); + return {apply_nm, D>(p1), p}; + } + +/** remove the rotational component of a, leaving only the poincare part */ +TD multivector despin(const multivector& a) { + auto p = a; + for(int i=(1<<(D::Dim-1)); i<(1<<(D::Dim)); i++) p[i] = typename D::Number(0); + auto p1 = chkmul,rotational,poincare>(a, conjugate(p)); + if(nm == nmInvariant) return p1 * pow(chkmul,rotational,units>(p, conjugate(p))[0], -0.5); + return apply_nm, D>(p1); + } + +TD std::string nzv(const mvector& a) { return "vector(" + nz(embed(a)) + ")"; } +TD std::string nzv(const matrix& a) { return ""; } + +template +typename D::Number sqnorm(multivector a) { + using N = typename D::Number; + auto res = chkmul>(a, conjugate(a))[0]; + if(res <= N(0) || isinf(res) || isnan(res)) res = N(1); + return res; + } + +TD typename D::Number sqnorm(mvector a) { + using N = typename D::Number; + N res(0); + for(int i=0; i multivector flatten(multivector a) { + using N = typename D::Number; + auto divby = a[0]; a[0] = N(1); + for(int i=1; i<(1< +multivector apply_nm(multivector a) { + if(nm == nmFlatten) return flatten(a); + if(nm == nmForced || nm == nmWeak) return a * pow(sqnorm(a), -0.5); + if(nm == nmBinary) { while(a[0] >= 2) { a = a / 2; } while(a[0] > 0 && a[0] < 0.5) { a = a * 2; } } + return a; + } + +TD mvector apply_nm(mvector a) { + if(nm == nmFlatten) { cbc[cbcDiv]--; return a / a[D::Dim-1]; } + if(nm == nmForced || nm == nmWeak) return a * pow(sqnorm(a), -0.5); + if(nm == nmBinary) { while(a[D::Dim-1] >= 2) { a = a / 2; } while(a[D::Dim-1] > 0 && a[D::Dim-1] < 0.5) { a = a * 2; } } + return a; + } + +/** get b which is a coordinate of a, but in normalized form. That is, if a is normalized simply return b, otherwise, multiply b appropriately */ + +template E get_normalized(multivector a, E b) { + if(nm != nmInvariant && nm != nmForced) return b * pow(sqnorm(a), -0.5); + return b; + } + +template E get_normalized(mvector a, E b) { + if(nm != nmInvariant && nm != nmForced) return b * pow(sqnorm(a), -0.5); + return b; + } + +} + diff --git a/devmods/reps/rep-halfplane.cpp b/devmods/reps/rep-halfplane.cpp new file mode 100644 index 00000000..d2efcd05 --- /dev/null +++ b/devmods/reps/rep-halfplane.cpp @@ -0,0 +1,182 @@ +/** representation based on the halfplane model; assumes Dim=3 */ + +namespace reps { + +template struct sl2 : public array { + + sl2(F a, F b, F c, F d) { self[0] = a; self[1] = b; self[2] = c; self[3] = d; } + + sl2 operator * (const sl2& sec) const { + return sl2( + self[0] * sec[0] + self[1] * sec[2], + self[0] * sec[1] + self[1] * sec[3], + self[2] * sec[0] + self[3] * sec[2], + self[2] * sec[1] + self[3] * sec[3] + ); + } + + std::string print() { + return hr::lalign(0, "[", self[0], ",", self[1], ";", self[2], ",", self[3], "]"); + } + + }; + +TD sl2 split_quaternion_to_sl2(const multivector& h) { + auto h3 = h[0], h2 = h[1 | 2], h1 = h[1 | 4], h0 = h[2 | 4]; + return sl2(h3 - h1, h2 + h0, -h2 + h0, h3 + h1); + } + +TD multivector sl2_to_split_quaternion(const sl2& e) { + auto h0 = (e[1] + e[2]) / 2; + auto h3 = (e[0] + e[3]) / 2; + auto h1 = (e[3] - e[0]) / 2; + auto h2 = (e[1] - e[2]) / 2; + auto res = zero_vector>(); + res[0] = h3; res[1 | 2] = h2; res[1 | 4] = h1; res[2 | 4] = h0; + return res; + } + +template using sl2c = sl2>; + +TD sl2c split_biquaternion_to_sl2c(const multivector& h) { + using cn = std::complex; + return sl2(cn(h[0]-h[9], h[15]-h[6]), cn(h[3]+h[10], -h[5]-h[12]), cn(h[10]-h[3], h[12]-h[5]), cn(h[0]+h[9], h[6]+h[15])); + } + +TD multivector sl2c_to_split_biquaternion(const sl2c& e) { + auto res = zero_vector>(); + res[0] = +(real(e[0]) + real(e[3])) / 2; + res[3] = +(real(e[1]) - real(e[2])) / 2; + res[5] = -(imag(e[1]) + imag(e[2])) / 2; + res[6] = +(imag(e[3]) - imag(e[0])) / 2; + res[9] = +(real(e[3]) - real(e[0])) / 2; + res[10] = +(real(e[1]) + real(e[2])) / 2; + res[12] = +(imag(e[2]) - imag(e[1])) / 2; + res[15] = +(imag(e[0]) + imag(e[3])) / 2; + return res; + } + +TD struct rep_halfplane { + + using data = D; + using N = typename D::Number; + using point = std::complex; + using isometry = sl2; + + static isometry cspin(int i, int j, N alpha) { + // return split_quaternion_to_sl2( rep_clifford::cspin(i, j, alpha) ); + if(i>j) std::swap(i, j), alpha = -alpha; alpha /= 2; + auto ca = cos(alpha), sa = sin(alpha); + return isometry(ca, -sa, sa, ca); + } + static isometry cspin90(int i, int j, N alpha) { + // return split_quaternion_to_sl2( rep_clifford::cspin(i, j, alpha) ); + auto ca = sqrt(N(2)), sa = sqrt(N(2)); + if(i>j) std::swap(i, j), sa = -sa; + return isometry(ca, -sa, sa, ca); + } + static isometry lorentz(int i, int j, N alpha) { + // return split_quaternion_to_sl2( rep_clifford::lorentz(i, j, alpha) ); + if(i>j) std::swap(i, j); alpha /= 2; + if(i == 0) return isometry(exp(-alpha), N(0), N(0), exp(alpha)); + if(i == 1) { + auto ca = cosh(alpha), sa = sinh(alpha); + return isometry(ca, sa, sa, ca); + } + throw hr::hr_exception("bad lorentz"); + } + static isometry id() { return isometry(N(1),N(0),N(0),N(1)); }; + static point center() { return point(N(0), N(1)); }; + static point apply(const isometry& T, const point& x) { + return (T[0] * x + T[1] * 1) / (T[2] * x + T[3] * 1); + }; + static isometry apply(const isometry& T, const isometry& U) { return T * U; }; + + static typename rep_clifford::point to_poincare(const point& x) { + auto a = real(x), b = imag(x); + + auto tmp = isometry(sqrt(b), a/sqrt(b), N(0), N(1)/sqrt(b)); + auto sq = sl2_to_split_quaternion(tmp); + + // sq[0] = (sqrt(b) + 1/sqrt(b)) / 2;; sq[1 | 2] = a/sqrt(b)/2; sq[1 | 4] = (1/sqrt(b) - sqrt(b)) / 2; sq[2 | 4] = a/sqrt(b)/2; + + sq = despin(sq); + return typename rep_clifford::point({{sq}}); + } + + static isometry inverse(isometry T) { return isometry(T[3], -T[1], -T[2], T[0]); } + static isometry push(const point& p) { return split_quaternion_to_sl2(to_poincare(p)[0]); } + + static N dist0(const point& x) { return rep_clifford::dist0(to_poincare(x)); } + static N angle(const point& x) { return rep_clifford::angle(to_poincare(x)); } + static N get_coord(const point& x, int i) { return rep_clifford::get_coord(to_poincare(x), i); } + + // imag may be very small and still important, so do not use the default complex print + static std::string print(const point& x) { return hr::lalign(0, "{real:", real(x), " imag:", imag(x), "}"); } + static std::string print(const isometry& x) { return x.print(); } + }; + +TD struct rep_halfspace { + + using data = D; + using N = typename D::Number; + struct point { std::complex xy; N z; }; + using isometry = sl2c; + + static isometry cspin(int i, int j, N alpha) { + return split_biquaternion_to_sl2c( rep_clifford::cspin(i, j, alpha) ); + } + static isometry cspin90(int i, int j) { + return split_biquaternion_to_sl2c( rep_clifford::cspin90(i, j) ); + } + static isometry lorentz(int i, int j, N alpha) { + return split_biquaternion_to_sl2c( rep_clifford::lorentz(i, j, alpha) ); + } + static isometry id() { return isometry(N(1),N(0),N(0),N(1)); } + static point center() { return point{ .xy = N(0), .z = N(1) }; } + static point apply(const isometry& T, const point& x) { + auto nom = T[0] * x.xy + T[1] * N(1); + auto nomz= T[0] * x.z; + auto den = T[2] * x.xy + T[3] * N(1); + auto denz= T[2] * x.z; + // D = den + denz * j + auto dnorm = std::norm(den) + std::norm(denz); + using std::conj; + // conj(D) = conj(den) - denz * j + // N / D = (nom + nomz * j) / (den + denz * j) = + // = (nom + nomz * j) * (conj(den) - denz * j) / dnorm + + // auto rxy = (nom * conj(den) - nomz * j * denz * j); + // auto rz*j = (-nom * denz * j + nomz * j * conj(den)) + + // apply the formula: j * a = conj(a) * j + + auto rxy = (nom * conj(den) + nomz * conj(denz)); + auto rz = (nomz * den - nom * denz); // todo only real part + // println(hlog, "imag of rz = ", imag(rz)); + return point { .xy = rxy / dnorm, .z = real(rz) / dnorm }; + }; + static isometry apply(const isometry& T, const isometry& U) { return T * U; }; + + static typename rep_clifford::point to_poincare(const point& x) { + auto tmp = isometry(sqrt(x.z), x.xy/sqrt(x.z), N(0), N(1)/sqrt(x.z)); + auto sq = sl2c_to_split_biquaternion(tmp); + sq = despin(sq); + return typename rep_clifford::point({{sq}}); + } + + static isometry inverse(isometry T) { return isometry(T[3], -T[1], -T[2], T[0]); } + static isometry push(const point& p) { return split_biquaternion_to_sl2c(to_poincare(p)[0]); } + + static N dist0(const point& x) { return rep_clifford::dist0(to_poincare(x)); } + static N angle(const point& x) { return rep_clifford::angle(to_poincare(x)); } + static N get_coord(const point& x, int i) { return rep_clifford::get_coord(to_poincare(x), i); } + + // imag may be very small and still important, so do not use the default complex print + static std::string print(const point& x) { return hr::lalign(0, "{x:", real(x.xy), " y:", imag(x.xy), " z:", x.z, "}"); } + static std::string print(const isometry& x) { return x.print(); } + }; + +template using rep_half = typename std::conditional, rep_halfspace>::type; + +} diff --git a/devmods/reps/rep-hr.cpp b/devmods/reps/rep-hr.cpp new file mode 100644 index 00000000..f90dd2ed --- /dev/null +++ b/devmods/reps/rep-hr.cpp @@ -0,0 +1,27 @@ +namespace reps { + +/* pull the HyperRogue representation; assumes HyperRogue geometry is set correctly, Number = ld, and Dim=3 or 4 */ +TD struct rep_hr { + + using data = D; + using N = typename D::Number; + using point = hr::hyperpoint; + using isometry = hr::transmatrix; + + static constexpr isometry cspin(int i, int j, N alpha) { return hr::cspin(i, j, ld(alpha)); } + static constexpr isometry cspin90(int i, int j) { return hr::cspin90(i, j); } + static constexpr isometry lorentz(int i, int j, N alpha) { return hr::lorentz(i, j, ld(alpha)); } + static isometry id() { return hr::Id; }; + static point center() { return D::Dim == 4 ? hr::C03 : hr::C02; }; + static point apply(const isometry& T, const point& x) { return T * x; }; + static isometry apply(const isometry& T, const isometry& U) { return T * U; }; + static ld dist0(const point& x) { return hdist0(x); } + static ld angle(const point& x) { return atan2(x[1], x[0]); } + static ld get_coord(const point& x, int i) { return x[i]; } + static isometry inverse(const isometry& T) { return iso_inverse(T); } + static isometry push(const point& p) { return rgpushxto0(p); } + + static std::string print(point p) { return hr::lalign(0, p); } + static std::string print(isometry p) { return hr::lalign(0, p); } + }; +} diff --git a/devmods/reps/rep-multi.cpp b/devmods/reps/rep-multi.cpp new file mode 100644 index 00000000..4b6cef30 --- /dev/null +++ b/devmods/reps/rep-multi.cpp @@ -0,0 +1,355 @@ +namespace reps { + +TD typename D::Number acos_auto(typename D::Number x) { + using N = typename D::Number; + if(hyperbolic) { + if(x < N(1)) return N(0); + return acosh(x); + } + if(sphere) { + if(x > N(1)) return N(0); + return acos(x); + } + throw hr::hr_exception("error"); + } + +/* use the linear representation, as in HyperRogue, but DO NOT apply nm, for comparison */ +TD struct rep_linear_nn { + + using data = D; + using point = mvector; + using isometry = matrix; + using N = typename D::Number; + + static constexpr isometry id() { + matrix result; + for(int i=0; i(D::Dim-1); } + static point apply(const isometry& T, const point& x) { return T * x; }; + static isometry apply(const isometry& T, const isometry& U) { return T * U; }; + static typename D::Number dist0(point x) { + return acos_auto (x[D::Dim-1]); + } + static typename D::Number angle(const point& x) { return atan2(x[1], x[0]); } + static typename D::Number get_coord(point x, int i) { return x[i]; } + + static isometry inverse(isometry T) { + for(int i=0; i get_column(matrix a, int id) { + mvector tmp; + for(int i=0; i a, mvector b) { + using N = typename D::Number; + N res(0); + for(int i=0; i& a, int id, mvector v) { + for(int i=0; i a) { return sqnorm(get_column(a, D::Dim-1)); } + +bool fix_matrices; + +TD matrix apply_nm(matrix a) { + using N = typename D::Number; + // normalize first + auto& lead = a[D::Dim-1][D::Dim-1]; + if(nm == nmFlatten) a = a / lead, cbc[cbcDiv]--; + if(nm == nmForced || nm == nmWeak) a = a * pow(sqnorm(a), -0.5); + if(nm == nmBinary) { + while(lead >= 2 && !isinf(lead)) { a = a / 2; } while(lead > 0 && lead < 0.5) { a = a * 2; } + } + // fixmatrix later + if(!fix_matrices) return a; + auto divby = (nm == nmBinary || nm == nmWeak || nm == nmCareless || nm == nmFlatten) ? sqnorm(a) : N(1); + for(int i=D::Dim-2; i>=0; i--) { + auto ro = get_column(a, i); + auto last = get_column(a, D::Dim-1); + ro = ro + last * inner(ro, last) / divby; + for(int j=i+1; j N(0)) ro = ro * (pow(in*in, -.5) * divby); + set_column(a, i, ro); + } + return a; + } + +/* use the linear representation, as in HyperRogue */ +TD struct rep_linear { + + using data = D; + using point = mvector; + using isometry = matrix; + using N = typename D::Number; + + static constexpr isometry cspin(int i, int j, typename D::Number angle) { + return apply_nm( rep_linear_nn::cspin(i, j, angle) ); + } + static constexpr isometry cspin90(int i, int j) { + return rep_linear_nn::cspin90(i, j); + } + static constexpr isometry lorentz(int i, int j, typename D::Number angle) { + return apply_nm( rep_linear_nn::lorentz(i, j, angle) ); + } + static isometry id() { return rep_linear_nn::id(); }; + static constexpr point center() { return unit_vector(D::Dim-1); } + static point apply(const isometry& T, const point& x) { return apply_nm(T * x); }; + static isometry apply(const isometry& T, const isometry& U) { return apply_nm(T * U); }; + static typename D::Number dist0(point x) { + return acos_auto (get_normalized(x, x[D::Dim-1])); + } + static typename D::Number angle(const point& x) { return atan2(x[1], x[0]); } + static typename D::Number get_coord(point x, int i) { + return get_normalized(x, x[i]); } + + static isometry inverse(isometry T) { + return rep_linear_nn::inverse(T); + } + static isometry push(const point& p) { + return apply_nm( rep_linear_nn::push(get_normalized(p, p)) ); + } + + static std::string print(point p) { return nzv(p); } + static std::string print(isometry p) { return nzv(p); } + }; + +/* use the linear representation of points and the multivector representation of isometries */ + +TD struct rep_mixed { + + using data = D; + using N = typename D::Number; + + using point = mvector; + using isometry = multivector; + + static isometry cspin(int i, int j, typename data::Number alpha, bool noflat = false) { +/* auto u = unit_vector> (0); + auto ui = unit_vector (i); + auto uj = unit_vector (j); + return u * cos(alpha/2) + multimul(embed(ui), embed(uj)) * sin(alpha/2); */ + auto res = zero_vector> (); + if(nm == nmFlatten && !noflat) { + res[0] = N(1); + res[(1< j ? 1 : -1); + return res; + } + res[0] = cos(alpha/2); + res[(1< j ? 1 : -1); + return res; + } + + static isometry cspin90(int i, int j, bool noflat = false) { + auto res = zero_vector> (); + if(nm == nmFlatten && !noflat) { + res[0] = N(1); + res[(1< j ? 1 : -1); + return res; + } + res[0] = sqrt(N(.5)); + res[(1< j ? 1 : -1); + return res; + } + + static isometry lorentz(int i, int j, typename data::Number alpha) { +/* // j must be time coordinate + auto u = unit_vector> (0); + auto ui = unit_vector (i); + auto uj = unit_vector (j); + return u * cosh(alpha/2) + multimul(embed(uj), embed(ui)) * sinh(alpha/2); */ + auto res = zero_vector> (); + if(nm == nmFlatten) { + res[0] = N(1); + res[(1<> (0); }; + static constexpr point center() { return unit_vector(D::Dim-1); } + static point apply(const isometry& T, const point& x) { + // return unembed(multimul(multimul(T, embed(x)), conjugate(T))); + return apply_nm(unembed(chkmul,flat_even,underling>(chkmul,flat_underling,odd>(T, embed(x)), conjugate(T)))); + }; + static isometry apply(const isometry& T, const isometry& U) { + auto res = apply_nm, D>(chkmul,flat_even,even>(T, U)); + return res; + } + + static isometry inverse(isometry T) { return conjugate(T); } + static isometry push(const point& p) { + auto pm = get_normalized(p, p); + pm[D::Dim-1] = pm[D::Dim-1] + N(1); + // since p was normalized, sqnorm of pm is 2 * pm[D::Dim-1] + pm = pm * pow(2 * pm[D::Dim-1], -0.5); + multivector v1 = embed(pm); + multivector v2 = unit_vector>(1<<(D::Dim-1)); + multivector v3 = chkmul,underling,poincare>(v1, v2); + v3 = apply_nm, D>(v3); + return v3; + } + + static typename D::Number dist0(point x) { return acos_auto (get_normalized(x, x[D::Dim-1])); } + static typename D::Number angle(const point& x) { return atan2(x[1], x[0]); } + static typename D::Number get_coord(point x, int i) { return get_normalized(x, x[i]); } + + static std::string print(point p) { return nzv(p); } + static std::string print(isometry p) { return nz(p); } + }; + +/* use the hyperboloid-Poincare representation of points and the multivector representation of isometries */ + +TD struct rep_clifford { + + using data = D; + using N = typename D::Number; + + using point = array< multivector, 1>; + using isometry = multivector; + + static isometry cspin(int i, int j, typename data::Number alpha) { return rep_mixed::cspin(i, j, alpha); } + static isometry cspin90(int i, int j) { return rep_mixed::cspin90(i, j); } + // j must be the neg coordinate! + static isometry lorentz(int i, int j, N alpha) { return rep_mixed::lorentz(i, j, alpha); } + static isometry id() { return rep_mixed::id(); } + static constexpr point center() { return point{{ id() }}; } + static point apply(const isometry& T, const point& x) { return point{{ despin(chkmul,poincare,even>(T, x[0])) }}; } + static isometry apply(const isometry& T, const isometry& U) { return apply_nm,D>( chkmul,even,even>(T, U) ); } + + static isometry inverse(isometry T) { return conjugate(T); } + static isometry push(const point& p) { return p[0]; } + + static typename D::Number dist0(const point& ax) { + return acos_auto(get_normalized, D, N>(ax[0], ax[0][0]))*2; + } + static constexpr int mvlast = 1<<(D::Dim-1); + static typename D::Number angle(const point& x) { + return atan2(x[0][2 | mvlast], x[0][1 | mvlast]); + } + static typename D::Number get_coord(const point& x, int i) { + auto x1 = multimul(multimul(x[0], unit_vector> (mvlast)), conjugate(x[0])); + auto x2 = unembed(x1); + return get_normalized(x2, x2[i]); + } + + static std::string print(point p) { return nz(p[0]); } + static std::string print(isometry p) { return nz(p); } + }; + +/* split isometries into the poincare and rotational part */ + +TD struct rep_gyro { + + using data = D; + using N = typename D::Number; + + using point = multivector; + using isometry = poincare_rotation; + + static isometry cspin(int i, int j, typename data::Number alpha) { return { rep_mixed::id(), rep_mixed::cspin(i, j, alpha, true) }; } + static isometry cspin90(int i, int j, typename data::Number alpha) { return { rep_mixed::id(), rep_mixed::cspin90(i, j, alpha, true) }; } + static isometry lorentz(int i, int j, typename data::Number alpha) { return {rep_mixed::lorentz(i, j, alpha), rep_mixed::id() }; } + static isometry id() { return { rep_mixed::id(), rep_mixed::id() }; } + static constexpr point center() { return rep_mixed::id(); } + static point apply(const isometry& T, const point& x) { return despin(chkmul,poincare,even>(T.first, chkmul,poincare,poincare>(T.second, x))); } + static isometry apply(const isometry& T, const isometry& U) { + auto R1 = apply_nm, poincare, poincare> (T.second, U.first); + auto R2 = apply_nm, poincare, even> (T.first, R1); + auto R3 = despin2(R2); + return { R3.first, apply_nm, rotational, rotational> (R3.second, U.second) }; + } + + static isometry inverse(isometry T) { return { conjugate(T.first), conjugate(T.second) }; } + static isometry push(const point& p) { return { p, rep_mixed::id() }; } + + static typename D::Number dist0(const point& ax) { + return acos_auto(get_normalized, D, N>(ax, ax[0]))*2; + } + static constexpr int mvlast = 1<<(D::Dim-1); + static typename D::Number angle(const point& x) { + return atan2(x[0][2 | mvlast], x[0][1 | mvlast]); + } + static typename D::Number get_coord(const point& x, int i) { + auto x1 = multimul(multimul(x[0], unit_vector> (mvlast)), conjugate(x[0])); + auto x2 = unembed(x1); + return get_normalized(x2, x2[i]); + } + + static std::string print(point p) { return nz(p[0]); } + static std::string print(isometry p) { return "["+nz(p.first)+","+nz(p.second)+"]"; } + }; + +} diff --git a/devmods/reps/rep-polar.cpp b/devmods/reps/rep-polar.cpp new file mode 100644 index 00000000..48d69d02 --- /dev/null +++ b/devmods/reps/rep-polar.cpp @@ -0,0 +1,261 @@ +bool polar_mod = true, polar_choose = true; + +namespace reps { + +template static void cyclefix(N& a) { + while(a > + get_deg(180)) a -= get_deg(360); + while(a < - get_deg(180)) a += get_deg(360); + } + +template static N cyclefix_on(N a) { cyclefix(a); return a; } + +/** the Taylor polynomial for 1-sqrt(1-y*y) */ +template N ssqrt(N y) { + return y*y/2 + y*y*y*y/8 + y*y*y*y*y*y*y/16 + y*y*y*y*y*y*y*y*y/128; + } + +TD struct rep_polar2 { + + using data = D; + using N = typename D::Number; + + struct point { N phi, r; }; + struct isometry { N psi, phi, r; }; // spin by psi first + + static isometry cspin(int i, int j, N alpha) { + if(i>j) std::swap(i, j), alpha = -alpha; + return isometry{.psi = -alpha, .phi = N(0), .r = N(0) }; + } + static isometry cspin90(int i, int j) { return cspin(i, j, get_deg(90)); } + static isometry lorentz(int i, int j, N alpha) { + if(i>j) std::swap(i, j); + if(i == 0) return isometry{.psi = N(0), .phi = N(0), .r = alpha}; + if(i == 1) return isometry{.psi = N(0), .phi = get_deg(90), .r = alpha}; + throw hr::hr_exception("bad lorentz"); + } + static isometry id() { return isometry{.psi = N(0), .phi = N(0), .r = N(0)}; }; + static point center() { return point{.phi = N(0), .r = N(0)}; }; + + static std::string print(isometry T) { + return hr::lalign(0, "{phi=", T.phi, " r=", T.r, " psi=", T.psi, "}"); + } + + static std::string print(point T) { + return hr::lalign(0, "{phi=", T.phi, " r=", T.r, "}"); + } + + static isometry apply(isometry T, isometry U, bool need_psi = true) { + if(T.r == 0) return isometry {.psi = T.psi+U.psi, .phi = T.psi+U.phi, .r = U.r}; + if(U.r == 0) return isometry {.psi = T.psi+U.psi, .phi = T.phi, .r = T.r}; + N alpha = U.phi + T.psi - T.phi; + if(polar_mod) cyclefix(alpha); + isometry res; + N y1 = sinh(U.r) * sin(alpha); + + auto ca = cos(alpha); + auto sa = sin(alpha); + N x1, x2; + + // choose the appropriate method + if(polar_choose && ca >= N(0.5)) { + N u = ca >= N(.999999) ? ssqrt(sa) : N(1) - ca; + res.r = cosh(T.r + U.r) - u * sinh(T.r) * sinh(U.r); + x1 = sinh(T.r + U.r) - u * cosh(T.r) * sinh(U.r); + if(need_psi) x2 = sinh(T.r + U.r) - u * cosh(U.r) * sinh(T.r); + } + else if(polar_choose && ca <= N(-0.5)) { + N u = ca <= N(-.999999) ? ssqrt(-sa) : ca + N(1); + res.r = cosh(T.r - U.r) + u * sinh(T.r) * sinh(U.r); + x1 = sinh(T.r - U.r) + u * cosh(T.r) * sinh(U.r); + if(need_psi) x2 = sinh(U.r - T.r) + u * cosh(U.r) * sinh(T.r); + } + else { + res.r = sinh(T.r) * sinh(U.r) * ca + cosh(T.r) * cosh(U.r); + x1 = cosh(T.r) * sinh(U.r) * ca + cosh(U.r) * sinh(T.r); + if(need_psi) x2 = cosh(U.r) * sinh(T.r) * ca + cosh(T.r) * sinh(U.r); + } + + if(res.r < N(1)) res.r = N(0); else res.r = acosh(res.r); + + N beta = (y1 || x1) ? atan2(y1, x1) : N(0); + res.phi = T.phi + beta; + if(polar_mod) cyclefix(res.phi); + + if(need_psi) { + N y2 = sinh(T.r) * sin(alpha); + N gamma = (y2 || x2) ? atan2(y2, x2) : N(0); + res.psi = T.psi + U.psi + beta + gamma - alpha; + if(polar_mod) cyclefix(res.psi); + } + + return res; + }; + + static point apply(const isometry& T, const point& x) { + isometry x1 = apply(T, push(x), false); + return point { .phi = x1.phi, .r = x1.r}; + }; + + static isometry inverse(isometry T) { return isometry{.psi = -T.psi, .phi = cyclefix_on(get_deg(180)+T.phi-T.psi), .r=T.r }; }; + static isometry push(const point& p) { return isometry{.psi = N(0), .phi = p.phi, .r = p.r}; } + + static N dist0(const point& x) { return x.r; } + static N angle(const point& x) { return x.phi; } + static N get_coord(const point& x, int i) { + if(i == 0) return cos(x.phi) * sinh(x.r); + if(i == 1) return sin(x.phi) * sinh(x.r); + if(i == 2) return cosh(x.r); + throw hr::hr_exception("bad get_coord"); + } + }; + +TD struct rep_high_polar { + + using data = D; + using N = typename D::Number; + + struct sphere_data { + using Number = N; + static constexpr int Dim = D::Dim-1; + static constexpr int Flipped = -1; + }; + + using subsphere = rep_linear_nn; + + struct point { typename subsphere::point phi; N r; }; + struct isometry { typename subsphere::isometry psi; typename subsphere::point phi; N r; }; + + static isometry cspin(int i, int j, N alpha) { + return isometry{.psi = subsphere::cspin(i, j, alpha), .phi = subsphere::center(), .r = N(0) }; + } + static isometry cspin90(int i, int j) { + return isometry{.psi = subsphere::cspin90(i, j), .phi = subsphere::center(), .r = N(0) }; + } + static isometry lorentz(int i, int j, N alpha) { + if(i>j) std::swap(i, j); + auto is = isometry{.psi = subsphere::id(), .phi = subsphere::center(), .r = alpha}; + is.phi[D::Dim-2] = N(0); + is.phi[i] = N(1); + return is; + } + static isometry id() { return isometry{.psi = subsphere::id(), .phi = subsphere::center(), .r = N(0)}; } + static point center() { return point{.phi = subsphere::center(), .r = N(0)}; }; + + static std::string print(isometry T) { + return hr::lalign(0, "{phi=", subsphere::print(T.phi), " r=", T.r, " psi=", hr::kz(T.psi.values), "}"); + } + + static std::string print(point T) { + return hr::lalign(0, "{phi=", subsphere::print(T.phi), " r=", T.r, "}"); + } + + static isometry apply(isometry T, isometry U, bool need_psi = true) { + auto apsi = need_psi ? T.psi * U.psi : subsphere::id(); + if(T.r == 0) return isometry {.psi = apsi, .phi = T.psi*U.phi, .r = U.r}; + if(U.r == 0) return isometry {.psi = apsi, .phi = T.phi, .r = T.r}; + auto aphi = T.psi * U.phi; + auto cos_alpha = inner(aphi, T.phi); + auto& ca = cos_alpha; + + isometry res; + + N x1, x2; + + auto orth = (aphi - T.phi * ca); + N sin_alpha; + if(ca > N(0.999999) || ca < N(0.999999)) + sin_alpha = pow(sqnorm(orth), .5); + else + sin_alpha = pow(N(1) - ca * ca, .5); + + if(sin_alpha == N(0)) { + if(ca >= N(1)) { + return isometry{.psi = apsi, .phi = T.phi, .r = T.r + U.r }; + } + if(ca <= N(-1)) { + if(T.r >= U.r) { + return isometry{.psi = apsi, .phi = T.phi, .r = T.r - U.r }; + } + else { + return isometry{.psi = apsi, .phi = T.phi*-1, .r = U.r - T.r }; + } + } + } + + orth = orth / sin_alpha; + N y1 = sinh(U.r) * sin_alpha; + + // choose the appropriate method + if(polar_choose && ca >= N(0.5)) { + N u = ca >= N(.999999) ? ssqrt(sin_alpha) : N(1) - ca; + res.r = cosh(T.r + U.r) - u * sinh(T.r) * sinh(U.r); + x1 = sinh(T.r + U.r) - u * cosh(T.r) * sinh(U.r); + if(need_psi) x2 = sinh(T.r + U.r) - u * cosh(U.r) * sinh(T.r); + } + else if(polar_choose && ca <= N(-0.5)) { + N u = ca <= N(-.999999) ? ssqrt(sin_alpha) : ca + N(1); // ca = u - 1 + res.r = cosh(T.r - U.r) + u * sinh(T.r) * sinh(U.r); + x1 = sinh(T.r - U.r) + u * cosh(T.r) * sinh(U.r); + if(need_psi) x2 = sinh(U.r - T.r) + u * cosh(U.r) * sinh(T.r); + } + else { + res.r = sinh(T.r) * sinh(U.r) * ca + cosh(T.r) * cosh(U.r); + x1 = cosh(T.r) * sinh(U.r) * ca + cosh(U.r) * sinh(T.r); + if(need_psi) x2 = cosh(U.r) * sinh(T.r) * ca + cosh(T.r) * sinh(U.r); + } + + if(res.r < N(1)) res.r = N(0); else res.r = acosh(res.r); + + auto h1 = pow(x1*x1+y1*y1, -0.5); + N cos_beta = x1*h1, sin_beta = y1*h1; + + res.phi = T.phi * cos_beta + orth * sin_beta; + + if(need_psi) { + N y2 = sinh(T.r) * sin_alpha; + auto h2 = pow(x2*x2+y2*y2, -0.5); + N cos_gamma = x2*h2, sin_gamma = y2*h2; + + // delta = beta + gamma - alpha + + auto cos_beta_gamma = cos_beta * cos_gamma - sin_beta * sin_gamma; + auto sin_beta_gamma = cos_beta * sin_gamma + sin_beta * cos_gamma; + + auto cos_delta = cos_beta_gamma * cos_alpha + sin_beta_gamma * sin_alpha; + auto sin_delta = sin_beta_gamma * cos_alpha - cos_beta_gamma * sin_alpha; + + auto phi1 = T.phi * cos_delta + orth * sin_delta; + auto orth1 = orth * cos_delta - T.phi * sin_delta; + + auto phi2 = phi1 - T.phi; + auto orth2 = orth1 - orth; + typename subsphere::isometry spinner = subsphere::id(); + + // Tv = v + * (phi'-phi) + * (orth'-orth) + + for(int i=0; i using rep_polar = typename std::conditional, rep_high_polar>::type; + +} diff --git a/devmods/reps/reps.cpp b/devmods/reps/reps.cpp new file mode 100644 index 00000000..26d23a91 --- /dev/null +++ b/devmods/reps/reps.cpp @@ -0,0 +1,81 @@ +#include + +#include "../../hyper.h" + +#define TD template +#undef sl2 + +namespace reps { + +using namespace boost::multiprecision; +using big = mpfr_float_50; +} + +namespace hr { + void print(hr::hstream& hs, ::reps::big b) { + std::stringstream ss; + ss << std::setprecision(10); + ss << b; string u; ss >> u; print(hs, u); + } +} + +namespace reps { +using std::array; +using std::vector; + +using hr::cell; +using hr::print; +using hr::hlog; +using hr::celldistance; +using hr::ld; +using hr::ginf; +using hr::geometry; +using hr::gcHyperbolic; +using hr::gcSphere; +using hr::C02; +using hr::C03; +using hr::qANYQ; + +template N get_deg(int deg); +template<> ld get_deg (int deg) { return M_PI*deg/180; } +template<> big get_deg (int deg) { return atan(big(1))*deg/45; } + +enum eNormalizeMode { + nmInvariant, // if the input was normalized, the output will be normalized too + nmForced, // normalize the output + nmWeak, // weakly normalize the output + nmCareless, // do not try to keep the output normalized + nmFlatten, // flatten the representation + nmBinary // try to avoid overflow + }; + +eNormalizeMode nm; + +} + +#include "counter.cpp" +#include "multivector.cpp" +#include "rep-hr.cpp" +#include "rep-multi.cpp" +#include "rep-halfplane.cpp" +#include "rep-polar.cpp" +#include "tests.cpp" + +namespace reps { + +// -- tests --- + +void test_systems() { + run_all_tests(); + fflush(stdout); + exit(1); + } + +void set_repgeo() { + if(test_dim == 3) { hr::set_geometry(hr::gNormal); hr::set_variation(hr::eVariation::pure); } + if(test_dim == 4) { hr::set_geometry(hr::gSpace435); } + } + +int a = hr::arg::add1("-test-reps", test_systems) + hr::arg::add1("-repgeo", set_repgeo); + +} diff --git a/devmods/reps/results.Md b/devmods/reps/results.Md new file mode 100644 index 00000000..ebe644e2 --- /dev/null +++ b/devmods/reps/results.Md @@ -0,0 +1,733 @@ +# What is it + +This is a study of numerical precision errors in various representations of 2D hyperbolic geometry. +It is generally the best to combine a representation with tiling; the tests take this into +account. + +# Representations studied + +The following representations are studied: + +* **linear**: points in the hyperboloid model; isometries as linear transformation matrices. + +* **mixed**: points in the hyperboloid model; isometries using Clifford algebras. (Clifford algebras + are a generalization of 'quaternions' commonly used in 3D graphics.) + +* **clifford**: points are also represented using Clifford algebras, that is, p is represented as + the isometry u such as u(C0) = p and u does not introduce extra rotations. + +* **halfplane (2D)**: points are represented using the half-plane model; isometries are represented using + SL(2,R). + +* **halfspace (3D)**: points are represented using the half-space model; isometries are represented using + SL(2,C). + +* **polar 2D**: points are represented using polar coordinates; isometries need one extra angle. + +* **general polar**: like polar 2D, but instead of angles, we use rotated unit vectors and rotation + matrices; this also makes it work in higher dimension. + +## Variations + +Both in linear and Clifford representations, there is the correct "normalized" representation; +if the normalized representation is multiplied by some factor x, most formulas still work, +and for those which do not, it is easy to compute x. This yields the following variations: + +* **invariant**: keep the invariant that the points and isometries are normalized + (that is: output is normalized under the assumption that the input is normalized) + +* **careless**: do not care about normalization + (advantages: some computations are avoided; possible to represent ultra-ideal points in + linear representations) + +* **forced**: normalize the output after every computation + (might be a good idea for points/isometries close to the center, but generally a bad + idea if they are far away -- in that case, the norm generally cannot be computed, but + distances and angles still tend to be correct in the invariant computations) + +* **weakly forced**: like forced, but do not normalize if the norm could not be computed + due to precision errors + +* **flatten**: instead of normal normalization, make the leading coordinate equal to 1. + The leading coordinate is the 'timelike' coordinate of linear representations + of points, and the 'unit' coordinate of Clifford representations. + (advantage: save memory: H2 represented only 2 coordinates instead of 3; + disadvantage: might not represent ultra-ideal points if they would be infinite) + +* **binary**: in careless, values may easily explode and cause underflow/overflow; avoid this + by making the leading coordinate in \[0.5, 2) range (by multiplying by powers of 2, which is + presumably fast) + +Furthermore: +* in linear, matrices can be **fixed** by replacing them by a correct orthogonal matrix close + to the current computation +* in (non-general) polar, forcing angles into [-pi,pi] may be needed to prevent explosion +* in **improved** polar, one of three variants of the cosine rule can be used, depending on the angle, + to improve the numerical precision; also even more precise computation to avoid numerical + precision errors for angles very close to 0 or pi +* in the Clifford representation, the **gyro** variant splits the isometries into + the translational part (which is flattened, making it equivalent to the Poincare disk model) + and the rotational part (for which 'invariant' is used). This fixes the problem + with full flattening where rotations by 180° are flattened to infinity. (AFAIK + Hyperbolica uses roughly this) + +## Observations + +* except linear, all the methods of representing isometries can only represent + orientation-preserving ones + +* Clifford isometries of H2 is essentially the same as SL(2,R) of halfplane -- it is + just the change of the basis + +* linear/Clifford representations are not that good at representing points close to the + boundary of the disk (invariant can somewhat tell the distance but flattened cannot); + halfplane is better here + +# Tests + +## test_loop_iso + +In this test, for each i, we construct a path in the tiling by always moving to a random +adjacent tile, until we get to a tile i afar; then, we return to the start (also randomly, +may stray further from the path). We compose all the relative tile isometries into T and see +if T(C0) = C0. The score is the first i for which it fails. + +Discussion: This makes rep_mixed worse than rep_lorentz. + +## test_loop_point + +Same as test_loop_iso but we apply the consecutive isometries to point right away. + +Discussion: This makes rep_mixed worse than rep_lorentz. + +## test_angledist + +For each i (skipping some), construct a path outwards in the tiling, compose isometries, +and see if the distance and angle to that tile have been computed correctly. + +Discussion: Invariant representations have no problem with this, even if the points obtained are beyond the precision otherwise. + +## test_similarity, test_dissimilarity, test_other + +For each i, compute the distance between two points in distance i from the starting point. +The angle between them is very small (test_similarity), close to 180° (test_dissimilarity), +close to 1° (test_other). + +Discussion: Similarity is obviously the most difficult. Halfplane is surprisingly powerful in all cases. + +## test_walk + +This is essentially walking in a straight line in HyperRogue. After some time, it can be often clearly observed that we +have 'deviated' from the original straight line. This test checks how long we can walk. + +We construct an isometry T representing a random direction. In each step, we compose this isometry with a translation (T := T * translate(1/16)). +Whenever the point T * C0 is closer to the center of another tile, we rebase to that new tile. + +For a test, we actually do this in parallel with two isometries T0 and T1, where T1 = T0 * translate(1/32). We count the number of steps +until the paths diverge. Numbers over 1000 are not actually that good, 1000+n means that, after n steps, the implementation no longer detects tile +changes. Numbers of 10000 signify that some even weirder problem happened. + +Discussion: Since the isometry matrices are always small (constrained to tiles), fixing definitely helps here. Without fixing, T stops +being an isometry (an effect visible in HyperRogue when fixing is disabled). + +## test_close + +Here we see whether small errors accumulate when moving close to the center. In test i, we move randomly until we reach distance i+1, +after which we return to the start (always reducing the distance). After each return to the start, we check if the representation is +still fine (if not, we restart with the original representation). The number given is the number of errors in 10000 steps. + +Discussion: Errors do not appear to accumulate when we simply move close to the start (or rather, they accumulate very slowly). + +## test_count + +This simply computes the number of numerical operations performed for every geometric operation. Numerical operations are categorized as: +* Addition/subtraction +* Multiplication (multiplication by constant not counted) +* Division (division by constant not counted) +* Functions: exp, log, (a)sin/cos/tan(h), sqrt, inverse sqrt + +Geometric operations are: +* spin: return rotation by given angle in given axes +* L0: return translation by given value in axis 0 +* L1: return translation by given value in axis 1 +* ip: apply isometry to point +* ii: compose isometries +* d0: compute the distance of point from 0 +* angle: compute the (2D) angle of point +* inverse: compute the inverse of an isometry +* push: convert a point into a translation + +# Implementation notes + +Note: the program currently assumes hyperbolic geometry (it was intended to support spherical +geometry but not everywhere the correct handling is implemented). + +# Results + +## Results on the {7,3} tiling + +``` + test_loop_iso + linear+F invariant: (17,16,17,16,17,17,17,17,16,17,17,17,17,16,16,17,17,16,17,17) + linear+F forced : (20,19,19,20,20,19,20,20,19,20,20,18,19,19,19,19,20,19,19,19) + linear+F weak : (21,19,20,20,20,19,24,20,19,25,24,18,19,19,19,19,20,19,19,19) + linear+F flatten : (19,19,21,20,20,20,20,20,19,19,19,21,19,20,19,19,19,19,19,20) + linear+F careless : (19,19,19,19,20,19,20,20,19,20,18,18,19,18,19,19,19,19,19,19) + linear+F binary : (19,19,19,19,20,19,20,20,19,20,18,18,19,18,19,19,19,19,19,19) + linear-F invariant: (17,17,17,18,18,18,21,17,16,17,19,18,17,18,23,17,18,17,19,17) + linear-F forced : (19,19,19,19,20,19,19,20,19,20,20,20,19,19,19,19,19,19,19,19) + linear-F weak : (19,19,20,19,20,19,19,20,19,21,20,20,21,19,19,19,20,19,19,19) + linear-F flatten : (20,19,20,20,21,19,20,18,19,19,20,21,19,21,20,19,19,19,19,20) + linear-F careless : (20,19,19,21,19,19,20,17,19,20,19,20,19,19,19,19,19,19,19,20) + linear-F binary : (20,19,19,21,19,19,20,17,19,20,19,20,19,19,19,19,19,19,19,20) + mixed invariant: (34,35,34,36,34,35,34,34,36,35,34,35,36,34,35,34,33,35,32,36) + mixed forced : (34,34,34,35,34,36,33,34,36,35,35,36,36,34,35,34,34,35,34,34) + mixed weak : (34,34,34,35,34,36,33,34,36,35,35,36,36,34,35,34,34,35,34,34) + mixed flatten : (20,5,18,13,13,13,4,13,9,29,25,9,36,22,19,30,5,35,14,2) + mixed careless : (34,34,34,35,34,34,35,34,36,35,35,36,36,34,35,34,35,35,34,34) + mixed binary : (34,34,34,35,34,34,35,34,36,35,35,36,36,34,35,34,35,35,34,34) + Clifford invariant: (32,35,34,36,34,35,34,34,36,33,34,35,36,34,35,34,33,34,32,36) + Clifford forced : (34,34,34,35,34,36,33,34,36,35,35,36,36,34,35,34,34,35,34,34) + Clifford weak : (34,34,34,35,34,36,33,34,36,35,35,36,36,34,35,34,34,35,34,34) + Clifford flatten : (34,34,34,35,34,36,35,34,36,35,35,36,37,34,35,34,36,35,34,34) + Clifford careless : (34,34,34,35,34,34,35,34,36,35,35,36,36,34,35,34,35,35,34,34) + Clifford binary : (34,34,34,35,34,34,35,34,36,35,35,36,36,34,35,34,35,35,34,34) + Clifford gyro : (34,34,34,35,34,36,35,34,36,35,35,36,36,34,35,34,35,35,34,34) + halfplane invariant: (34,34,34,35,34,36,35,34,36,35,35,36,36,34,35,34,35,35,34,34) + polar basic : (34,34,34,35,34,34,33,34,36,35,35,24,36,34,35,34,35,35,34,34) + polar improved : (34,34,34,34,33,36,35,33,36,35,35,36,35,34,35,34,34,35,34,35) + polar F/F : (34,36,35,35,33,34,33,33,36,35,35,36,36,34,35,34,35,35,34,34) + polar F/T : (34,34,34,34,33,36,35,33,36,35,35,36,35,34,35,34,34,35,34,35) + polar T/F : (34,34,35,36,34,36,35,34,35,35,35,36,35,34,35,34,35,35,34,36) + polar T/T : (34,36,35,34,34,36,35,34,36,35,35,36,35,34,35,34,34,35,34,35) + test_loop_point + linear+F invariant: (17,16,17,16,17,17,17,17,16,17,17,17,17,16,16,17,17,16,17,17) + linear+F forced : (20,19,19,20,20,19,20,20,19,20,20,18,19,19,19,19,20,19,19,19) + linear+F weak : (21,19,20,20,20,19,24,20,19,25,24,18,19,19,19,19,20,19,19,19) + linear+F flatten : (19,19,21,20,20,20,20,20,19,19,19,21,19,20,19,19,19,19,19,20) + linear+F careless : (19,19,19,19,20,19,20,20,19,20,18,18,19,18,19,19,19,19,19,19) + linear+F binary : (19,19,19,19,20,19,20,20,19,20,18,18,19,18,19,19,19,19,19,19) + linear-F invariant: (17,17,17,18,18,18,21,17,16,17,19,18,17,18,23,17,18,17,19,17) + linear-F forced : (19,19,19,19,20,19,19,20,19,20,20,20,19,19,19,19,19,19,19,19) + linear-F weak : (19,19,20,19,20,19,19,20,19,21,20,20,21,19,19,19,20,19,19,19) + linear-F flatten : (20,19,20,20,21,19,20,18,19,19,20,21,19,21,20,19,19,19,19,20) + linear-F careless : (20,19,19,21,19,19,20,17,19,20,19,20,19,19,19,19,19,19,19,20) + linear-F binary : (20,19,19,21,19,19,20,17,19,20,19,20,19,19,19,19,19,19,19,20) + mixed invariant: (18,17,17,19,20,18,19,19,17,17,18,17,17,16,16,18,25,17,17,19) + mixed forced : (19,19,19,19,19,19,19,17,19,20,19,20,19,19,19,20,19,19,20,20) + mixed weak : (19,23,19,19,19,19,21,17,19,23,19,20,24,20,19,20,22,19,23,20) + mixed flatten : (20,19,19,20,19,19,19,18,20,20,19,18,18,19,18,20,19,19,19,19) + mixed careless : (19,19,19,19,19,20,19,18,19,19,19,21,19,20,18,19,20,19,19,20) + mixed binary : (19,19,19,19,19,20,19,18,19,19,19,21,19,20,18,19,20,19,19,20) + Clifford invariant: (32,34,31,34,33,36,31,35,32,33,32,36,32,34,33,36,33,32,34,34) + Clifford forced : (34,34,34,35,34,36,35,34,35,35,35,36,35,34,35,34,35,35,34,34) + Clifford weak : (34,34,34,35,34,36,35,34,35,35,35,36,35,34,35,34,35,35,34,34) + Clifford flatten : (34,34,34,35,34,36,35,34,36,35,34,36,36,34,35,34,35,35,34,34) + Clifford careless : (3,1,3,2,2,2,2,2,2,3,2,2,3,2,2,3,2,2,1,3) + Clifford binary : (34,34,34,35,34,34,35,34,36,35,35,36,35,34,35,34,35,35,34,34) + Clifford gyro : (34,36,34,35,34,36,35,34,36,35,35,36,36,34,35,34,35,35,34,34) + halfplane invariant: (34,36,34,35,34,36,35,34,36,35,35,36,36,34,35,34,35,35,34,34) + polar basic : (34,34,34,35,34,34,33,34,36,35,35,24,36,34,35,34,35,35,34,34) + polar improved : (34,34,34,34,33,36,35,33,36,35,35,36,35,34,35,34,34,35,34,35) + polar F/F : (34,36,35,35,33,34,33,33,36,35,35,36,36,34,35,34,35,35,34,34) + polar F/T : (34,34,34,34,33,36,35,33,36,35,35,36,35,34,35,34,34,35,34,35) + polar T/F : (34,34,35,36,34,36,35,34,35,35,35,36,35,34,35,34,35,35,34,36) + polar T/T : (34,36,35,34,34,36,35,34,36,35,35,36,35,34,35,34,34,35,34,35) + test_angledist + linear+F invariant: (767,767,767) + linear+F forced : (21,21,21) + linear+F weak : (21,21,21) + linear+F flatten : (21,21,21) + linear+F careless : (21,21,21) + linear+F binary : (21,21,21) + linear-F invariant: (767,767,767) + linear-F forced : (21,21,21) + linear-F weak : (21,21,21) + linear-F flatten : (21,21,21) + linear-F careless : (21,21,21) + linear-F binary : (21,21,21) + mixed invariant: (767,767,767) + mixed forced : (21,21,21) + mixed weak : (21,21,21) + mixed flatten : (21,21,21) + mixed careless : (21,21,21) + mixed binary : (21,21,21) + Clifford invariant: (767,767,767) + Clifford forced : (39,39,47) + Clifford weak : (39,39,39) + Clifford flatten : (39,39,39) + Clifford careless : (2,3,3) + Clifford binary : (39,47,39) + Clifford gyro : (39,47,39) + halfplane invariant: (767,767,767) + polar basic : (443,443,443) + polar improved : (443,443,443) + polar F/F : (767,767,767) + polar F/T : (767,767,767) + polar T/F : (767,767,767) + polar T/T : (767,767,767) + test_similarity + linear+F invariant: (18,17,17,18,17,18,18,17,18,17,18,18,17,17,17,18,17,17,17,17) + linear+F forced : (19,18,18,18,19,18,18,19,19,18,18,18,19,18,19,19,18,19,19,18) + linear+F weak : (19,18,18,18,19,18,18,19,19,18,18,18,19,18,19,19,18,19,19,18) + linear+F flatten : (18,19,18,18,18,19,19,19,19,18,18,18,19,18,19,19,19,19,18,19) + linear+F careless : (19,18,19,18,19,18,18,19,19,18,18,18,19,19,19,19,18,18,19,18) + linear+F binary : (19,18,19,18,19,18,18,19,19,18,18,18,19,19,19,19,18,18,19,18) + linear-F invariant: (18,18,19,18,18,18,18,18,18,18,18,19,19,19,18,18,19,18,18,18) + linear-F forced : (18,19,19,19,19,19,18,19,18,19,19,19,19,19,18,19,19,19,19,19) + linear-F weak : (18,19,19,19,19,19,18,19,18,19,19,19,19,19,18,19,19,19,19,19) + linear-F flatten : (19,18,19,19,19,19,19,19,18,19,19,19,19,19,18,19,19,19,19,19) + linear-F careless : (18,19,19,19,19,19,19,19,18,19,19,19,19,19,19,19,19,19,19,19) + linear-F binary : (18,19,19,19,19,19,19,19,18,19,19,19,19,19,19,19,19,19,19,19) + mixed invariant: (18,19,18,19,18,18,18,18,18,19,19,19,18,19,18,18,18,19,19,18) + mixed forced : (19,19,19,19,19,18,19,19,18,18,18,19,20,19,19,20,19,19,19,19) + mixed weak : (19,19,19,19,19,18,19,19,18,18,19,19,21,20,19,19,19,20,19,20) + mixed flatten : (19,19,19,19,19,19,19,19,20,18,19,19,19,19,20,19,19,19,19,19) + mixed careless : (19,20,19,19,19,18,19,19,19,19,19,20,20,19,19,19,19,19,19,19) + mixed binary : (19,20,19,19,19,18,19,19,19,19,19,20,20,19,19,19,19,19,19,19) + Clifford invariant: (33,33,33,32,32,34,34,33,32,34,33,33,34,33,33,33,33,34,33,33) + Clifford forced : (36,35,35,35,35,35,37,35,36,36,36,35,35,35,35,37,36,36,35,35) + Clifford weak : (36,35,35,35,35,35,37,35,36,36,36,35,35,35,35,37,36,36,35,35) + Clifford flatten : (35,37,36,35,36,35,35,38,36,37,36,35,37,35,37,38,36,38,35,36) + Clifford careless : (37,35,36,35,35,36,35,35,37,35,38,36,38,35,37,36,35,35,35,37) + Clifford binary : (37,35,36,35,35,36,35,35,37,35,38,36,38,35,37,36,35,35,35,37) + Clifford gyro : (37,35,37,37,37,36,35,36,36,36,37,36,35,36,37,37,36,37,37,35) + halfplane invariant: (35,36,36,36,36,35,35,37,37,37,37,36,35,36,38,37,36,36,37,35) + polar basic : (19,18,18,18,19,18,18,18,18,18,18,18,18,19,18,18,18,18,18,18) + polar improved : (37,37,37,37,37,38,38,37,37,36,37,37,38,39,38,37,37,39,37,38) + polar F/F : (19,19,19,19,19,19,19,19,19,19,19,19,19,19,19,19,19,19,19,19) + polar F/T : (35,35,35,35,35,35,35,35,35,35,35,35,35,35,35,35,35,35,35,35) + polar T/F : (19,19,19,19,19,19,19,19,19,19,19,19,19,19,19,19,19,19,19,19) + polar T/T : (35,35,35,35,35,36,35,35,35,36,36,35,35,35,35,35,35,35,35,36) + test_dissimilarity + linear+F invariant: (125,124,147,123,134,130,126,128,123,125,130,130,127,125,125,131,123,124,125,127) + linear+F forced : (7,6,7,7,7,7,8,7,6,7,7,7,7,7,7,7,7,7,7,7) + linear+F weak : (7,7,7,7,9,13,8,7,9,7,9,7,7,7,11,9,7,8,7,7) + linear+F flatten : (7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7) + linear+F careless : (7,7,7,7,7,7,7,7,6,7,8,7,7,7,7,7,7,7,7,7) + linear+F binary : (7,7,7,7,7,7,7,7,6,7,8,7,7,7,7,7,7,7,7,7) + linear-F invariant: (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + linear-F forced : (10,10,10,10,10,10,10,10,10,9,9,10,10,10,10,10,10,9,10,10) + linear-F weak : (9,15,13,12,11,11,10,9,16,9,10,9,13,10,12,11,10,9,10,9) + linear-F flatten : (9,10,10,10,10,10,9,9,10,10,10,10,10,10,9,10,10,10,10,10) + linear-F careless : (10,10,10,9,10,9,10,10,10,10,9,10,10,9,9,10,9,10,10,9) + linear-F binary : (10,10,10,9,10,9,10,10,10,10,9,10,10,9,9,10,9,10,10,9) + mixed invariant: (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + mixed forced : (10,9,10,9,10,10,10,9,10,10,10,10,10,10,10,10,10,9,10,10) + mixed weak : (10,9,12,9,9,10,10,11,10,9,12,11,10,9,11,11,11,9,11,10) + mixed flatten : (10,9,10,10,9,10,10,10,10,10,10,10,10,10,9,10,10,10,10,10) + mixed careless : (10,10,10,10,10,10,10,10,10,10,10,10,9,9,9,10,10,10,10,9) + mixed binary : (10,10,10,10,10,10,10,10,10,10,10,10,9,9,9,10,10,10,10,9) + Clifford invariant: (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + Clifford forced : (18,19,18,18,19,18,19,18,19,19,18,19,18,17,19,19,19,19,19,19) + Clifford weak : (18,18,18,18,19,18,18,18,19,18,18,19,18,18,18,18,19,18,18,19) + Clifford flatten : (18,19,18,19,19,18,18,18,19,19,18,18,19,19,18,18,18,19,19,18) + Clifford careless : (19,18,19,18,18,19,19,18,18,18,18,18,19,18,19,19,19,18,19,18) + Clifford binary : (19,18,19,18,18,19,19,18,18,18,18,18,19,18,19,19,19,18,19,18) + Clifford gyro : (18,18,18,19,19,18,19,20,18,18,18,18,18,19,19,19,18,19,18,18) + halfplane invariant: (35,35,35,35,35,35,35,35,34,36,37,34,35,36,35,35,35,35,36,35) + polar basic : (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + polar improved : (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + polar F/F : (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + polar F/T : (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + polar T/F : (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + polar T/T : (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + test_other + linear+F invariant: (97,97,101,107,97,97,126,107,101,101,112,112,99,101,112,131,107,97,97,94) + linear+F forced : (10,10,10,10,10,10,10,10,11,10,10,10,10,11,10,11,10,10,10,10) + linear+F weak : (10,11,10,11,12,13,11,11,13,14,10,10,10,11,11,11,12,10,13,10) + linear+F flatten : (10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,11,10) + linear+F careless : (10,11,10,10,10,10,10,10,10,10,10,10,10,11,10,10,10,10,10,10) + linear+F binary : (10,11,10,10,10,10,10,10,10,10,10,10,10,11,10,10,10,10,10,10) + linear-F invariant: (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + linear-F forced : (12,12,12,13,13,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12) + linear-F weak : (12,12,12,18,13,12,12,14,14,14,12,12,12,12,13,12,12,12,12,12) + linear-F flatten : (12,12,12,12,12,12,12,12,12,12,13,12,12,12,12,12,12,12,12,12) + linear-F careless : (12,12,13,12,12,12,13,12,13,13,12,12,12,12,12,12,12,12,12,12) + linear-F binary : (12,12,13,12,12,12,13,12,13,13,12,12,12,12,12,12,12,12,12,12) + mixed invariant: (359,360,359,359,359,359,359,359,359,359,359,359,359,359,360,359,359,359,359,359) + mixed forced : (14,13,14,14,13,13,13,14,13,14,14,14,14,13,14,13,13,14,13,13) + mixed weak : (13,13,14,14,13,13,13,14,13,18,13,16,14,13,15,14,13,14,13,13) + mixed flatten : (13,14,14,14,14,14,13,14,14,14,13,14,13,14,14,13,13,14,14,13) + mixed careless : (14,13,14,14,13,14,14,14,14,14,14,14,13,14,13,14,13,14,13,13) + mixed binary : (14,13,14,14,13,14,14,14,14,14,14,14,13,14,13,14,13,14,13,13) + Clifford invariant: (361,361,361,361,361,361,361,361,361,361,361,361,361,361,361,361,361,361,361,361) + Clifford forced : (21,21,23,21,21,21,22,21,23,21,22,21,21,21,21,22,21,21,21,22) + Clifford weak : (21,21,23,21,21,21,22,21,23,21,21,21,21,21,21,22,21,21,21,22) + Clifford flatten : (21,21,22,21,22,22,21,23,21,21,22,22,22,22,22,21,21,22,22,22) + Clifford careless : (21,22,22,21,21,22,22,23,21,21,21,22,22,21,21,21,22,22,22,23) + Clifford binary : (21,22,22,21,21,22,22,23,21,21,21,22,22,21,21,21,22,22,22,23) + Clifford gyro : (23,23,23,23,24,23,23,24,23,24,23,23,23,23,23,23,23,23,23,23) + halfplane invariant: (35,35,35,35,36,38,37,35,36,36,37,36,35,36,35,36,36,38,38,35) + polar basic : (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + polar improved : (360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360) + polar F/F : (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + polar F/T : (360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360) + polar T/F : (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + polar T/T : (360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360) + test_walk + linear+F invariant: (658,606,622,603,606,608,667,620,615,616,616,619,631,609,613,614,632,635,592,612) + linear+F forced : (632,617,606,617,608,607,607,628,627,632,626,627,596,652,623,639,615,617,615,616) + linear+F weak : (613,626,622,636,606,605,640,594,609,615,592,639,625,607,613,600,635,620,622,604) + linear+F flatten : (605,618,608,625,677,649,612,607,614,649,621,609,602,599,615,645,609,597,597,638) + linear+F careless : (617,615,620,608,603,606,654,597,612,598,626,624,601,613,616,609,602,627,601,619) + linear+F binary : (617,615,620,608,603,606,654,597,612,598,626,624,601,613,616,609,602,627,601,619) + linear-F invariant: (309,298,341,329,349,304,302,301,292,314,305,310,322,307,314,1315,1312,295,10000,10000) + linear-F forced : (1294,1290,1318,1308,1297,1287,1297,1295,1348,1304,1301,1309,1303,1318,292,1299,1294,1291,1288,1298) + linear-F weak : (298,295,301,289,299,289,287,301,297,294,305,297,316,286,292,291,301,293,303,299) + linear-F flatten : (297,1324,1299,1293,1315,1301,1313,1305,1305,306,1282,1310,1307,1309,1292,1306,1296,1315,1313,1289) + linear-F careless : (1305,1310,1299,1298,1308,307,1318,1296,1293,1300,302,1306,1298,1287,1299,1315,1312,290,1290,1298) + linear-F binary : (1305,1310,1299,1298,1308,307,1318,1296,1293,1300,302,1306,1298,1287,1299,1315,1312,290,1290,1298) + mixed invariant: (622,634,587,620,607,619,609,612,599,646,607,623,616,590,615,592,637,636,659,613) + mixed forced : (599,666,588,608,607,602,630,671,601,602,618,630,618,601,599,599,614,601,596,617) + mixed weak : (599,666,588,608,607,602,630,671,601,602,618,630,618,601,599,599,614,601,596,617) + mixed flatten : (611,616,607,605,610,603,595,605,593,614,617,593,602,642,610,616,625,593,636,617) + mixed careless : (622,634,587,620,607,619,609,612,599,646,607,623,616,590,615,592,637,636,659,613) + mixed binary : (622,634,587,620,607,619,609,612,599,646,607,623,616,590,615,592,637,636,659,613) + Clifford invariant: (622,634,587,620,607,619,609,612,599,646,607,623,616,590,615,592,637,636,659,613) + Clifford forced : (599,666,588,608,607,602,630,671,601,602,618,630,618,601,599,599,614,601,596,617) + Clifford weak : (599,666,588,608,607,602,630,671,601,602,618,630,618,601,599,599,614,601,596,617) + Clifford flatten : (611,616,607,605,610,603,595,605,593,614,617,593,602,642,610,616,625,593,636,617) + Clifford careless : (622,634,587,620,607,619,609,612,599,646,607,623,616,590,615,592,637,636,659,613) + Clifford binary : (622,634,587,620,607,619,609,612,599,646,607,623,616,590,615,592,637,636,659,613) + Clifford gyro : (600,586,602,621,625,621,625,603,593,630,634,600,586,597,600,609,601,592,617,615) + halfplane invariant: (600,586,602,621,625,621,625,603,593,630,634,600,586,597,600,609,601,592,617,615) + polar basic : (1055,67,1078,66,72,1050,66,1073,70,1052,66,1070,66,67,71,73,1064,1070,67,66) + polar improved : (71,1068,1073,1059,67,55,67,1071,65,1052,1067,1078,67,63,69,1067,57,66,69,1059) + polar F/F : (605,566,605,563,566,583,565,591,578,616,591,568,601,569,584,559,621,579,589,601) + polar F/T : (573,596,633,569,581,590,565,588,590,581,600,614,597,571,595,619,576,573,582,631) + polar T/F : (594,662,662,598,591,686,590,610,593,592,588,588,600,581,598,572,618,578,589,588) + polar T/T : (583,594,601,586,570,601,594,579,585,581,582,614,649,614,674,639,588,580,587,588) + test_close + linear+F invariant: (0,0,0,0,0,0,0,0,0,0,0,0,0,0,3,18,78,130,126,118) + linear+F forced : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,11,70,110,118) + linear+F weak : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,11,70,110,118) + linear+F flatten : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,10,55,99,116) + linear+F careless : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,15,67,107,113) + linear+F binary : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,15,67,107,113) + linear-F invariant: (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,27,89,117,118) + linear-F forced : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,19,73,117,117) + linear-F weak : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,19,73,117,116) + linear-F flatten : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,15,49,103,115) + linear-F careless : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,17,62,112,115) + linear-F binary : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,17,62,112,115) + mixed invariant: (0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,2,36,101,117,115) + mixed forced : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,21,72,114,117) + mixed weak : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,21,72,114,115) + mixed flatten : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,20,75,115,116) + mixed careless : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,17,73,113,112) + mixed binary : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,17,73,113,112) + Clifford invariant: (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + Clifford forced : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + Clifford weak : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + Clifford flatten : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + Clifford careless : (1666,1698,1241,859,666,545,447,378,339,298,262,245,244,207,196,175,170,168,157,144) + Clifford binary : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + Clifford gyro : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + halfplane invariant: (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + polar basic : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + polar improved : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + polar F/F : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + polar F/T : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + polar T/F : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + polar T/T : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + test_count + linear+F invariant: (spin(24A 34M 9D 4F) L0(24A 34M 9D 4F) L1(24A 34M 9D 4F) ip(9A 9M) ii(51A 51M 9D) d0(1F) angle(1F) inverse() push(29A 42M 10D 2F)) + linear+F forced : (spin(27A 46M 9D 5F) L0(27A 46M 9D 5F) L1(27A 46M 9D 5F) ip(12A 15M 1F) ii(54A 73M 9D 3F) d0(1F) angle(1F) inverse() push(32A 54M 10D 3F)) + linear+F weak : (spin(30A 49M 9D 5F) L0(30A 49M 9D 5F) L1(30A 49M 9D 5F) ip(12A 15M 1F) ii(57A 76M 9D 3F) d0(3A 4M 2F) angle(1F) inverse() push(38A 63M 10D 4F)) + linear+F flatten : (spin(27A 37M 17D 4F) L0(27A 37M 17D 4F) L1(27A 37M 17D 4F) ip(9A 9M 2D) ii(54A 54M 17D) d0(3A 4M 2F) angle(1F) inverse() push(35A 51M 18D 3F)) + linear+F careless : (spin(27A 37M 9D 4F) L0(27A 37M 9D 4F) L1(27A 37M 9D 4F) ip(9A 9M) ii(54A 64M 9D 2F) d0(3A 4M 2F) angle(1F) inverse() push(35A 51M 10D 3F)) + linear+F binary : (spin(27A 37M 9D 4F) L0(27A 37M 9D 4F) L1(27A 37M 9D 4F) ip(9A 9M) ii(54A 64M 9D 2F) d0(3A 4M 2F) angle(1F) inverse() push(35A 51M 10D 3F)) + linear-F invariant: (spin(2F) L0(2F) L1(2F) ip(9A 9M) ii(27A 27M) d0(1F) angle(1F) inverse() push(5A 8M 1D)) + linear-F forced : (spin(3A 12M 3F) L0(3A 12M 3F) L1(3A 12M 3F) ip(12A 15M 1F) ii(30A 39M 1F) d0(1F) angle(1F) inverse() push(8A 20M 1D 1F)) + linear-F weak : (spin(3A 12M 3F) L0(3A 12M 3F) L1(3A 12M 3F) ip(12A 15M 1F) ii(30A 39M 1F) d0(3A 4M 2F) angle(1F) inverse() push(11A 26M 1D 2F)) + linear-F flatten : (spin(8D 2F) L0(8D 2F) L1(8D 2F) ip(9A 9M 2D) ii(27A 27M 8D) d0(3A 4M 2F) angle(1F) inverse() push(8A 14M 9D 1F)) + linear-F careless : (spin(2F) L0(2F) L1(2F) ip(9A 9M) ii(27A 27M) d0(3A 4M 2F) angle(1F) inverse() push(8A 14M 1D 1F)) + linear-F binary : (spin(2F) L0(2F) L1(2F) ip(9A 9M) ii(27A 27M) d0(3A 4M 2F) angle(1F) inverse() push(8A 14M 1D 1F)) + mixed invariant: (spin(2F) L0(2F) L1(2F) ip(17A 24M) ii(12A 16M) d0(1F) angle(1F) inverse() push(5A 7M)) + mixed forced : (spin(2F) L0(2F) L1(2F) ip(20A 30M 1F) ii(15A 28M 1F) d0(1F) angle(1F) inverse() push(7A 18M 1F)) + mixed weak : (spin(2F) L0(2F) L1(2F) ip(20A 30M 1F) ii(15A 28M 1F) d0(3A 4M 2F) angle(1F) inverse() push(10A 24M 2F)) + mixed flatten : (spin(1F) L0(1F) L1(1F) ip(17A 15M 2D) ii(12A 12M) d0(3A 4M 2F) angle(1F) inverse() push(8A 15M 1F)) + mixed careless : (spin(2F) L0(2F) L1(2F) ip(17A 24M) ii(12A 16M) d0(3A 4M 2F) angle(1F) inverse() push(8A 13M 1F)) + mixed binary : (spin(2F) L0(2F) L1(2F) ip(17A 24M) ii(12A 16M) d0(3A 4M 2F) angle(1F) inverse() push(8A 13M 1F)) + Clifford invariant: (spin(2F) L0(2F) L1(2F) ip(12A 28M 1F) ii(12A 16M) d0(1F) angle(1F) inverse() push()) + Clifford forced : (spin(2F) L0(2F) L1(2F) ip(13A 29M 1F) ii(15A 28M 1F) d0(1F) angle(1F) inverse() push()) + Clifford weak : (spin(2F) L0(2F) L1(2F) ip(13A 29M 1F) ii(15A 28M 1F) d0(2A 4M 2F) angle(1F) inverse() push()) + Clifford flatten : (spin(1F) L0(1F) L1(1F) ip(11A 20M) ii(12A 19M) d0(2A 4M 2F) angle(1F) inverse() push()) + Clifford careless : (spin(2F) L0(2F) L1(2F) ip(11A 18M) ii(12A 16M) d0(2A 4M 2F) angle(1F) inverse() push()) + Clifford binary : (spin(2F) L0(2F) L1(2F) ip(11A 18M) ii(12A 16M) d0(2A 4M 2F) angle(1F) inverse() push()) + Clifford gyro : (spin(2F) L0(2F) L1(2F) ip(5A 10M 2D) ii(4A 8M) d0(9A 12M 2D 5F) angle(7A 8M 2D 4F) inverse() push(11A 8M 2D 3F)) + halfplane invariant: (spin(2F) L0(2F) L1(2F) ip(5A 10M 2D) ii(4A 8M) d0(8A 16M 2D 5F) angle(8A 16M 2D 5F) inverse() push(12A 16M 2D 4F)) + polar basic : (spin(2F) L0() L1() ip(15A 25M 2D 12F) ii(53A 73M 2D 17F) d0() angle(1F) inverse(4A 4M) push()) + polar improved : (spin(2F) L0() L1() ip(20A 41M 2D 10F) ii(59A 88M 2D 15F) d0() angle(1F) inverse(4A 4M) push()) + polar F/F : (spin() L0() L1() ip(5A 7M 14F) ii(10A 11M 21F) d0() angle() inverse(14A) push()) + polar F/T : (spin() L0() L1() ip(5A 7M 14F) ii(14A 8M 18F) d0() angle() inverse(14A) push()) + polar T/F : (spin() L0() L1() ip(5A 7M 14F) ii(11A 11M 21F) d0() angle() inverse(2A) push()) + polar T/T : (spin() L0() L1() ip(5A 7M 14F) ii(15A 8M 18F) d0() angle() inverse(2A) push()) +``` + +## Results on the {4,3,5} honeycomb + +``` + test_loop_iso + linear+F invariant: (32,32,34,34,41,39,60,43,32,31,38,36,36,36,41,50,33,37,33,53) + linear+F forced : (22,22,24,25,25,23,25,25,25,26,23,26,24,26,24,25,24,24,24,26) + linear+F weak : (24,25,26,25,28,26,25,25,27,25,23,27,27,27,27,26,24,24,27,26) + linear+F flatten : (22,22,24,25,25,23,27,26,25,26,25,26,24,26,26,26,24,24,24,26) + linear+F careless : (26,25,26,25,25,26,27,25,26,26,23,26,27,27,26,26,25,24,27,26) + linear+F binary : (26,25,26,25,25,26,27,25,26,26,23,26,27,27,26,26,25,24,27,26) + linear-F invariant: (35,72,27,37,38,43,27,25,29,44,24,47,26,27,28,29,42,35,29,29) + linear-F forced : (26,25,26,25,25,26,25,25,26,27,25,26,27,27,26,27,24,25,24,26) + linear-F weak : (31,25,30,25,26,29,25,25,26,27,25,26,29,27,27,27,24,31,24,26) + linear-F flatten : (22,22,26,25,25,23,25,26,25,26,23,24,24,26,24,27,24,24,24,26) + linear-F careless : (22,22,24,25,26,23,25,25,25,26,23,24,27,26,24,26,24,24,27,26) + linear-F binary : (22,22,24,25,26,23,25,25,25,26,23,24,27,26,24,26,24,24,27,26) + mixed invariant: (49,47,44,47,42,44,47,45,46,47,45,49,46,50,45,49,42,48,49,45) + mixed forced : (51,50,49,49,47,47,46,47,48,47,48,49,50,50,49,47,52,47,48,49) + mixed weak : (51,50,49,49,47,47,46,47,48,47,48,49,50,50,49,47,52,47,48,49) + mixed flatten : (1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1) + mixed careless : (51,50,46,50,47,47,48,47,48,47,48,49,50,50,49,51,52,51,48,49) + mixed binary : (51,50,46,50,47,47,48,47,48,47,48,49,50,50,49,51,52,51,48,49) + Clifford invariant: (46,47,44,44,42,44,47,45,45,46,45,50,46,50,45,39,42,44,49,45) + Clifford forced : (51,50,49,49,47,47,46,47,48,47,48,49,50,50,49,47,52,47,48,49) + Clifford weak : (51,50,49,49,47,47,46,47,48,47,48,49,50,50,49,47,52,47,48,49) + Clifford flatten : (55,78,999,58,57,81,55,79,98,999,58,51,66,52,99,51,58,54,88,66) + Clifford careless : (51,50,46,50,47,47,48,47,48,47,48,49,50,50,49,51,52,51,48,49) + Clifford binary : (51,50,46,50,47,47,48,47,48,47,48,49,50,50,49,51,52,51,48,49) + Clifford gyro : (51,50,49,50,50,47,48,47,51,47,50,49,50,50,49,51,52,54,50,51) + halfplane invariant: (51,50,49,50,47,47,48,47,48,47,50,49,50,50,49,51,52,51,48,51) + polar basic : (12,50,46,49,47,47,48,47,48,3,46,49,50,50,47,7,52,3,48,49) + polar improved : (49,46,49,50,50,47,48,47,51,47,50,49,50,50,49,51,48,51,48,49) + test_loop_point + linear+F invariant: (999,999,999,999,999,999,999,999,999,999,999,999,999,999,999,999,999,999,999,999) + linear+F forced : (22,22,24,25,25,23,25,25,25,26,23,26,24,26,24,25,24,24,24,26) + linear+F weak : (24,25,26,25,28,26,25,25,27,25,23,27,27,27,27,26,24,24,27,26) + linear+F flatten : (22,22,24,25,25,23,27,26,25,26,25,27,24,26,26,26,24,24,24,26) + linear+F careless : (26,25,26,25,25,26,27,25,26,26,23,26,27,27,26,26,25,24,27,26) + linear+F binary : (26,25,26,25,25,26,27,25,26,26,23,26,27,27,26,26,25,24,27,26) + linear-F invariant: (35,72,27,37,38,43,27,25,29,44,24,47,26,27,28,29,42,35,29,29) + linear-F forced : (26,25,26,25,25,26,25,25,26,27,25,26,27,27,26,27,24,25,24,26) + linear-F weak : (31,25,30,25,26,29,25,25,26,27,25,26,29,27,27,27,24,31,24,26) + linear-F flatten : (22,22,26,25,25,23,25,26,25,26,23,24,24,26,24,27,24,24,24,26) + linear-F careless : (22,22,24,25,26,23,25,25,25,26,23,24,27,26,24,26,24,24,27,26) + linear-F binary : (22,22,24,25,26,23,25,25,25,26,23,24,27,26,24,26,24,24,27,26) + mixed invariant: (24,22,22,25,23,25,26,22,24,22,23,25,26,27,24,23,25,22,23,24) + mixed forced : (26,25,26,25,25,23,25,25,26,27,23,26,27,26,26,26,25,24,24,26) + mixed weak : (27,25,30,25,25,23,25,25,26,27,23,26,27,26,26,28,25,24,24,28) + mixed flatten : (1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1) + mixed careless : (26,25,24,25,25,26,27,25,26,26,26,26,27,27,26,26,24,25,24,26) + mixed binary : (26,25,24,25,25,26,27,25,26,26,26,26,27,27,26,26,24,25,24,26) + Clifford invariant: (46,45,48,44,42,44,45,48,47,48,46,41,45,49,48,39,50,48,51,45) + Clifford forced : (49,50,46,49,47,47,46,45,48,47,45,49,50,50,49,51,52,47,48,49) + Clifford weak : (49,50,46,49,47,47,46,45,48,47,45,49,50,50,49,51,52,47,48,49) + Clifford flatten : (50,48,46,49,47,46,48,47,48,47,45,47,50,50,49,47,42,47,48,49) + Clifford careless : (3,3,1,2,2,2,2,2,3,3,3,2,3,2,1,2,3,3,1,3) + Clifford binary : (51,48,49,49,47,47,48,47,48,47,45,49,50,50,49,47,52,47,48,49) + Clifford gyro : (51,50,46,50,50,47,48,48,51,47,50,49,50,50,49,51,52,51,48,49) + halfplane invariant: (49,50,46,49,47,47,48,47,51,47,50,49,50,50,49,51,52,51,48,49) + polar basic : (12,50,46,49,47,47,48,47,48,3,46,49,50,50,47,7,52,3,48,49) + polar improved : (49,46,49,50,50,47,48,47,51,47,50,49,50,50,49,51,48,51,48,49) + test_angledist + linear+F invariant: (999,999,999) + linear+F forced : (26,26,32) + linear+F weak : (26,32,32) + linear+F flatten : (26,26,32) + linear+F careless : (26,26,32) + linear+F binary : (26,26,32) + linear-F invariant: (999,999,999) + linear-F forced : (26,26,32) + linear-F weak : (32,26,32) + linear-F flatten : (26,26,32) + linear-F careless : (26,32,32) + linear-F binary : (26,32,32) + mixed invariant: (999,999,999) + mixed forced : (26,26,32) + mixed weak : (26,26,32) + mixed flatten : (2,2,1) + mixed careless : (26,26,32) + mixed binary : (26,26,32) + Clifford invariant: (999,999,999) + Clifford forced : (57,57,47) + Clifford weak : (69,57,47) + Clifford flatten : (57,57,47) + Clifford careless : (5,4,4) + Clifford binary : (57,57,47) + Clifford gyro : (57,69,57) + halfplane invariant: (999,999,999) + polar basic : (532,532,5) + polar improved : (532,532,532) + test_similarity + linear+F invariant: (17,16,16,16,16,16,17,16,17,16,17,17,16,16,16,17,17,16,16,17) + linear+F forced : (17,17,17,18,19,18,18,18,18,17,18,18,17,18,17,18,18,17,17,17) + linear+F weak : (17,17,17,18,18,18,18,18,17,17,18,18,17,18,17,18,17,17,17,17) + linear+F flatten : (17,17,17,18,18,18,18,18,17,17,18,18,17,18,17,18,18,17,17,17) + linear+F careless : (17,17,17,18,19,18,17,18,17,17,18,18,17,18,17,17,18,17,17,17) + linear+F binary : (17,17,17,18,19,18,17,18,17,17,18,18,17,18,17,17,18,17,17,17) + linear-F invariant: (17,18,18,17,17,17,17,18,17,18,17,18,17,18,17,17,17,17,17,17) + linear-F forced : (18,19,18,18,18,18,18,18,18,18,18,19,19,18,19,19,18,19,18,18) + linear-F weak : (18,19,18,18,18,17,18,18,18,18,19,19,19,18,19,19,18,19,18,18) + linear-F flatten : (19,18,19,18,18,18,18,20,18,19,18,18,20,19,19,19,18,19,18,18) + linear-F careless : (19,18,18,18,19,18,18,20,19,18,18,19,19,19,18,18,19,18,18,18) + linear-F binary : (19,18,18,18,19,18,18,20,19,18,18,19,19,19,18,18,19,18,18,18) + mixed invariant: (18,19,19,19,19,19,18,19,19,18,18,19,19,19,19,19,19,18,19,18) + mixed forced : (19,19,19,18,18,19,19,18,19,19,19,19,19,19,19,20,19,19,19,19) + mixed weak : (19,19,19,18,18,19,19,18,18,19,19,19,19,20,19,18,19,18,18,19) + mixed flatten : (19,19,18,19,19,18,19,20,18,19,19,20,19,19,20,18,19,19,19,20) + mixed careless : (19,19,19,19,19,20,20,19,19,19,19,19,19,19,19,18,19,19,19,19) + mixed binary : (19,19,19,19,19,20,20,19,19,19,19,19,19,19,19,18,19,19,19,19) + Clifford invariant: (34,33,34,35,33,34,33,32,34,34,33,34,34,34,34,35,33,34,33,34) + Clifford forced : (35,36,35,37,36,35,35,36,36,36,35,37,35,36,36,36,36,36,36,35) + Clifford weak : (35,36,35,37,36,35,35,36,36,36,35,37,35,36,36,36,36,36,36,35) + Clifford flatten : (35,36,35,35,36,36,37,37,36,36,36,36,36,36,38,35,37,36,37,36) + Clifford careless : (36,36,38,36,37,35,35,37,36,36,37,35,37,35,35,35,35,37,36,38) + Clifford binary : (36,36,38,36,37,35,35,37,36,36,37,35,37,35,35,35,35,37,36,38) + Clifford gyro : (36,37,38,36,35,37,35,36,38,36,36,36,36,36,36,37,35,34,35,36) + halfplane invariant: (36,35,36,36,37,38,38,36,37,36,37,36,35,36,35,37,35,36,36,35) + polar basic : (17,18,18,17,17,17,17,18,18,18,17,17,17,18,17,17,17,17,17,17) + polar improved : (34,37,35,36,37,36,34,36,35,35,35,34,36,35,36,36,35,35,35,35) + test_dissimilarity + linear+F invariant: (63,63,64,75,63,63,64,64,64,64,63,63,73,66,63,63,76,74,64,64) + linear+F forced : (6,7,7,7,7,7,7,7,7,7,7,7,6,7,7,7,7,7,7,6) + linear+F weak : (7,9,8,8,7,8,8,7,13,8,7,8,6,7,7,8,8,7,8,7) + linear+F flatten : (8,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,6,7,7,7) + linear+F careless : (7,7,6,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,6) + linear+F binary : (7,7,6,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,6) + linear-F invariant: (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + linear-F forced : (10,10,10,10,10,10,10,10,10,10,10,9,9,10,9,10,10,9,10,9) + linear-F weak : (12,9,9,14,10,9,9,11,11,11,11,13,9,12,10,13,9,10,15,10) + linear-F flatten : (10,10,10,9,9,9,10,10,10,10,10,10,10,9,9,10,10,10,10,10) + linear-F careless : (10,9,9,10,9,10,9,10,9,10,10,10,10,10,9,9,9,10,10,9) + linear-F binary : (10,9,9,10,9,10,9,10,9,10,10,10,10,10,9,9,9,10,10,9) + mixed invariant: (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + mixed forced : (10,10,10,10,10,10,10,10,9,10,10,10,10,10,9,9,9,10,9,10) + mixed weak : (10,16,11,15,9,12,10,9,10,14,10,10,10,11,11,11,11,10,11,11) + mixed flatten : (9,10,9,10,10,10,9,10,9,10,10,10,10,9,10,9,10,9,9,9) + mixed careless : (10,9,9,10,10,10,10,9,10,10,9,10,10,10,10,10,9,9,9,10) + mixed binary : (10,9,9,10,10,10,10,9,10,10,9,10,10,10,10,10,9,9,9,10) + Clifford invariant: (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + Clifford forced : (18,17,17,19,20,17,18,18,19,19,18,18,17,19,18,18,19,18,19,18) + Clifford weak : (19,17,17,18,19,17,18,18,17,18,18,18,17,18,18,17,18,18,19,19) + Clifford flatten : (19,18,18,18,19,18,18,18,18,18,18,19,19,18,18,19,18,19,18,18) + Clifford careless : (18,18,18,18,18,18,18,19,19,18,18,18,18,18,18,19,18,19,19,18) + Clifford binary : (18,18,18,18,18,18,18,19,19,18,18,18,18,18,18,19,18,19,19,18) + Clifford gyro : (18,17,18,19,19,19,18,19,18,19,18,18,19,19,19,19,18,18,19,17) + halfplane invariant: (34,35,36,35,35,35,35,34,35,35,36,35,35,37,36,37,35,35,34,35) + polar basic : (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + polar improved : (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + test_other + linear+F invariant: (73,50,65,68,65,50,68,50,64,67,65,63,50,66,63,50,68,50,64,67) + linear+F forced : (10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10) + linear+F weak : (10,10,11,10,11,17,10,12,13,12,11,11,10,10,10,11,11,16,10,14) + linear+F flatten : (10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10) + linear+F careless : (10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10) + linear+F binary : (10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10) + linear-F invariant: (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + linear-F forced : (12,12,13,12,12,12,12,12,12,12,12,13,12,12,12,12,12,12,12,12) + linear-F weak : (12,12,12,12,13,12,14,13,12,15,12,14,16,12,15,12,12,16,12,12) + linear-F flatten : (12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12) + linear-F careless : (13,12,13,12,12,12,12,12,12,12,12,12,12,12,12,12,13,12,12,12) + linear-F binary : (13,12,13,12,12,12,12,12,12,12,12,12,12,12,12,12,13,12,12,12) + mixed invariant: (359,359,359,359,359,359,359,359,359,359,359,359,359,359,359,359,359,359,359,359) + mixed forced : (14,14,14,13,14,13,14,13,14,14,13,14,14,13,13,14,13,14,14,13) + mixed weak : (15,13,14,13,14,14,16,13,16,14,14,14,13,13,13,14,13,13,14,15) + mixed flatten : (14,13,13,14,13,14,14,14,13,13,14,14,14,14,13,14,14,13,14,14) + mixed careless : (14,14,13,13,13,13,13,13,14,14,13,14,14,13,13,14,13,13,14,13) + mixed binary : (14,14,13,13,13,13,13,13,14,14,13,14,14,13,13,14,13,13,14,13) + Clifford invariant: (361,361,361,361,361,361,361,361,361,361,361,361,361,361,361,361,361,361,361,361) + Clifford forced : (22,21,21,21,22,21,22,22,22,22,22,21,21,22,22,22,22,21,21,22) + Clifford weak : (22,21,21,21,22,21,22,22,22,22,22,21,21,22,22,22,22,21,21,22) + Clifford flatten : (22,22,22,22,22,22,22,22,22,22,22,21,21,22,22,22,22,22,22,22) + Clifford careless : (22,21,22,22,22,23,21,21,21,22,22,21,22,21,22,22,22,22,22,22) + Clifford binary : (22,21,22,22,22,23,21,21,21,22,22,21,22,21,22,22,22,22,22,22) + Clifford gyro : (23,24,23,23,23,24,23,23,23,23,23,23,23,23,23,23,23,24,23,23) + halfplane invariant: (35,35,35,35,35,35,35,35,35,36,37,35,35,35,36,36,35,35,35,35) + polar basic : (356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356,356) + polar improved : (360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360,360) + test_walk + linear+F invariant: (614,576,600,601,581,579,609,600,596,565,596,610,583,588,595,606,591,592,579,614) + linear+F forced : (598,653,583,603,584,597,581,582,618,568,608,621,586,596,636,604,584,619,585,586) + linear+F weak : (606,604,616,593,604,592,583,611,598,576,597,579,603,584,575,591,584,587,600,591) + linear+F flatten : (619,599,580,615,608,615,600,593,601,577,615,603,579,590,635,589,589,586,560,589) + linear+F careless : (620,596,575,593,602,591,576,587,603,574,596,604,589,615,599,623,586,620,612,602) + linear+F binary : (620,596,575,593,602,591,576,587,603,574,596,604,589,615,599,623,586,620,612,602) + linear-F invariant: (293,291,289,304,1295,307,285,302,290,297,1311,289,320,292,305,289,285,2132,2137,284) + linear-F forced : (1292,288,1303,1298,301,1288,1285,1290,297,288,304,298,293,299,1300,1291,1294,1276,1281,1288) + linear-F weak : (288,288,284,300,301,291,286,294,297,288,304,288,293,299,300,286,291,10000,281,293) + linear-F flatten : (1292,1297,1303,1311,1292,1285,1296,1297,10000,1299,1321,298,1311,1288,293,1315,1287,1280,1283,1286) + linear-F careless : (1288,1292,1288,303,1295,1295,1285,1303,291,1290,1310,289,1321,1290,1298,1289,1284,1276,1282,1287) + linear-F binary : (1288,1292,1288,303,1295,1295,1285,1303,291,1290,1310,289,1321,1290,1298,1289,1284,1276,1282,1287) + mixed invariant: (592,582,592,572,580,598,589,597,594,567,595,600,605,591,602,611,596,580,589,587) + mixed forced : (592,587,597,579,582,627,592,603,603,565,596,576,589,591,596,606,592,594,591,592) + mixed weak : (592,587,597,579,582,627,592,603,603,565,596,576,589,591,596,606,592,594,591,592) + mixed flatten : (28,363,52,303,296,336,292,316,288,287,308,298,298,298,310,314,293,371,353,324) + mixed careless : (592,582,592,572,580,598,589,597,594,567,595,600,605,591,602,611,596,580,589,587) + mixed binary : (592,582,592,572,580,598,589,597,594,567,595,600,605,591,602,611,596,580,589,587) + Clifford invariant: (592,582,592,572,580,598,589,597,594,567,595,600,605,591,602,611,596,580,589,587) + Clifford forced : (592,587,597,579,582,627,592,603,603,565,596,576,589,591,596,606,592,594,591,592) + Clifford weak : (592,587,597,579,582,627,592,603,603,565,596,576,589,591,596,606,592,594,591,592) + Clifford flatten : (614,604,595,598,570,595,589,587,599,573,599,585,603,602,634,602,596,600,584,596) + Clifford careless : (592,582,592,572,580,598,589,597,594,567,595,600,605,591,602,611,596,580,589,587) + Clifford binary : (592,582,592,572,580,598,589,597,594,567,595,600,605,591,602,611,596,580,589,587) + Clifford gyro : (601,581,606,584,595,604,609,589,582,596,601,601,579,601,619,589,588,583,588,587) + halfplane invariant: (576,614,594,596,580,604,594,600,610,566,593,584,603,597,606,608,588,579,575,590) + polar basic : (25,1047,69,1045,1045,1033,59,1065,1041,67,1064,56,1044,1034,70,68,55,53,53,55) + polar improved : (55,58,70,1045,1045,1057,56,1065,59,68,1065,55,1044,1034,1070,66,1056,53,52,58) + test_close + linear+F invariant: (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,4,22) + linear+F forced : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + linear+F weak : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + linear+F flatten : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + linear+F careless : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + linear+F binary : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + linear-F invariant: (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,8) + linear-F forced : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + linear-F weak : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + linear-F flatten : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + linear-F careless : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1) + linear-F binary : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1) + mixed invariant: (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1) + mixed forced : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1) + mixed weak : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1) + mixed flatten : (5000,2716,1760,1229,986,774,667,584,500,460,430,365,344,299,286,283,257,237,231,234) + mixed careless : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + mixed binary : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + Clifford invariant: (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + Clifford forced : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + Clifford weak : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + Clifford flatten : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + Clifford careless : (1666,1539,1427,1047,817,653,547,483,418,381,350,308,283,260,247,230,212,202,188,183) + Clifford binary : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + Clifford gyro : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + halfplane invariant: (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + polar basic : (0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + polar improved : (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) + test_count + linear+F invariant: (spin(60A 78M 24D 5F) L0(60A 78M 24D 5F) L1(60A 78M 24D 5F) ip(16A 16M) ii(124A 142M 24D 3F) d0(1F) angle(1F) inverse() push(70A 96M 25D 3F)) + linear+F forced : (spin(64A 98M 24D 6F) L0(64A 98M 24D 6F) L1(64A 98M 24D 6F) ip(20A 24M 1F) ii(128A 162M 24D 4F) d0(1F) angle(1F) inverse() push(74A 116M 25D 4F)) + linear+F weak : (spin(68A 102M 24D 6F) L0(68A 102M 24D 6F) L1(68A 102M 24D 6F) ip(20A 24M 1F) ii(132A 166M 24D 4F) d0(4A 5M 2F) angle(1F) inverse() push(82A 128M 25D 5F)) + linear+F flatten : (spin(64A 82M 39D 5F) L0(64A 82M 39D 5F) L1(64A 82M 39D 5F) ip(16A 16M 3D) ii(128A 146M 39D 3F) d0(4A 5M 2F) angle(1F) inverse() push(78A 108M 40D 4F)) + linear+F careless : (spin(64A 82M 24D 5F) L0(64A 82M 24D 5F) L1(64A 82M 24D 5F) ip(16A 16M) ii(128A 146M 24D 3F) d0(4A 5M 2F) angle(1F) inverse() push(78A 108M 25D 4F)) + linear+F binary : (spin(64A 82M 24D 5F) L0(64A 82M 24D 5F) L1(64A 82M 24D 5F) ip(16A 16M) ii(128A 146M 24D 3F) d0(4A 5M 2F) angle(1F) inverse() push(78A 108M 25D 4F)) + linear-F invariant: (spin(2F) L0(2F) L1(2F) ip(16A 16M) ii(64A 64M) d0(1F) angle(1F) inverse() push(10A 18M 1D)) + linear-F forced : (spin(4A 20M 3F) L0(4A 20M 3F) L1(4A 20M 3F) ip(20A 24M 1F) ii(68A 84M 1F) d0(1F) angle(1F) inverse() push(14A 38M 1D 1F)) + linear-F weak : (spin(4A 20M 3F) L0(4A 20M 3F) L1(4A 20M 3F) ip(20A 24M 1F) ii(68A 84M 1F) d0(4A 5M 2F) angle(1F) inverse() push(18A 46M 1D 2F)) + linear-F flatten : (spin(15D 2F) L0(15D 2F) L1(15D 2F) ip(16A 16M 3D) ii(64A 64M 15D) d0(4A 5M 2F) angle(1F) inverse() push(14A 26M 16D 1F)) + linear-F careless : (spin(2F) L0(2F) L1(2F) ip(16A 16M) ii(64A 64M) d0(4A 5M 2F) angle(1F) inverse() push(14A 26M 1D 1F)) + linear-F binary : (spin(2F) L0(2F) L1(2F) ip(16A 16M) ii(64A 64M) d0(4A 5M 2F) angle(1F) inverse() push(14A 26M 1D 1F)) + mixed invariant: (spin(2F) L0(2F) L1(2F) ip(52A 64M) ii(56A 64M) d0(1F) angle(1F) inverse() push(7A 10M)) + mixed forced : (spin(2F) L0(2F) L1(2F) ip(56A 72M 1F) ii(63A 88M 1F) d0(1F) angle(1F) inverse() push(10A 30M 1F)) + mixed weak : (spin(2F) L0(2F) L1(2F) ip(56A 72M 1F) ii(63A 88M 1F) d0(4A 5M 2F) angle(1F) inverse() push(14A 38M 2F)) + mixed flatten : (spin(1F) L0(1F) L1(1F) ip(52A 49M 3D) ii(56A 56M) d0(4A 5M 2F) angle(1F) inverse() push(11A 21M 1F)) + mixed careless : (spin(2F) L0(2F) L1(2F) ip(52A 64M) ii(56A 64M) d0(4A 5M 2F) angle(1F) inverse() push(11A 18M 1F)) + mixed binary : (spin(2F) L0(2F) L1(2F) ip(52A 64M) ii(56A 64M) d0(4A 5M 2F) angle(1F) inverse() push(11A 18M 1F)) + Clifford invariant: (spin(2F) L0(2F) L1(2F) ip(39A 68M 1F) ii(56A 64M) d0(1F) angle(1F) inverse() push()) + Clifford forced : (spin(2F) L0(2F) L1(2F) ip(39A 68M 1F) ii(63A 88M 1F) d0(1F) angle(1F) inverse() push()) + Clifford weak : (spin(2F) L0(2F) L1(2F) ip(39A 68M 1F) ii(63A 88M 1F) d0(3A 5M 2F) angle(1F) inverse() push()) + Clifford flatten : (spin(1F) L0(1F) L1(1F) ip(36A 51M) ii(56A 71M) d0(3A 5M 2F) angle(1F) inverse() push()) + Clifford careless : (spin(2F) L0(2F) L1(2F) ip(36A 48M) ii(56A 64M) d0(3A 5M 2F) angle(1F) inverse() push()) + Clifford binary : (spin(2F) L0(2F) L1(2F) ip(36A 48M) ii(56A 64M) d0(3A 5M 2F) angle(1F) inverse() push()) + Clifford gyro : (spin(8A 1F) L0(8A 1F) L1(8A 1F) ip(23A 38M 1D) ii(24A 32M) d0(23A 26M 1D 5F) angle(20A 21M 1D 4F) inverse() push(28A 21M 1D 3F)) + halfplane invariant: (spin(8A 2F) L0(8A 2F) L1(8A 2F) ip(23A 38M 1D) ii(24A 32M) d0(23A 38M 1D 5F) angle(23A 38M 1D 5F) inverse() push(31A 38M 1D 4F)) + polar basic : (spin(2F) L0() L1() ip(24A 35M 3D 12F) ii(114A 135M 3D 17F) d0() angle(1F) inverse(9A 9M) push()) + polar improved : (spin(2F) L0() L1() ip(29A 51M 3D 10F) ii(120A 150M 3D 15F) d0() angle(1F) inverse(9A 9M) push()) +``` \ No newline at end of file diff --git a/devmods/reps/tests.cpp b/devmods/reps/tests.cpp new file mode 100644 index 00000000..ea512a6b --- /dev/null +++ b/devmods/reps/tests.cpp @@ -0,0 +1,538 @@ +namespace reps { + +constexpr int test_dim = 3; +constexpr bool in_hyperbolic = true; + +int edges, valence; + +void prepare_tests() { + hr::start_game(); + if(MDIM != test_dim) throw hr::hr_exception("fix your dimension"); + if(!(in_hyperbolic ? hyperbolic : sphere)) throw hr::hr_exception("fix your geometry"); + if(hr::variation != hr::eVariation::pure) throw hr::hr_exception("fix your variation"); + if(quotient) throw hr::hr_exception("fix your quotient"); + if(test_dim == 4) { + if(cginf.tiling_name != "{4,3,5}") throw hr::hr_exception("only {4,3,5} implemented in 3D"); + edges = 4; + valence = 5; + } + else { + edges = hr::cwt.at->type; + bool ok; + valence = hr::get_valence(hr::cwt.at, 1, ok); + } + } + +struct data { + using Number = hr::ld; + static constexpr int Dim = test_dim; + static constexpr int Flipped = in_hyperbolic ? test_dim-1 : -1; + }; + +struct countdata { + using Number = countfloat; + static constexpr int Dim = data::Dim; + static constexpr int Flipped = data::Flipped; + }; + +struct bigdata { + using Number = big; + static constexpr int Dim = data::Dim; + static constexpr int Flipped = data::Flipped; + }; + +using good = rep_linear_nn; + +int debug; // 0 -- never, 1 -- only errors, 2 -- always + +vector randomwalk(std::mt19937& gen, cell *from, int dist) { + vector res = { from }; + while(celldistance(from, res.back()) < dist) { + int i = gen() % res.back()->type; + res.push_back(res.back()->cmove(i)); + } + return res; + } + +template N rand01(std::mt19937& gen) { return N(((gen() & HRANDMAX) | 1) / (HRANDMAX+1.)); } + +vector random_return(std::mt19937& gen, cell *from, cell *to, ld peq, ld pbad) { + vector res = { from }; + ld d = celldistance(to, from); + while(to != res.back()) { + int i = gen() % res.back()->type; + cell *r1 = res.back()->cmove(i); + ld d1 = celldistance(to, r1); + bool ok = d1 < d ? true : d1 == d ? rand01(gen) < peq : rand01(gen) < pbad; + if(ok) { res.push_back(r1); d = d1; } + } + return res; + } + +vector vrev(vector a) { reverse(a.begin(), a.end()); return a; } +vector vcon(vector a, vector b) { for(auto bi: b) a.push_back(bi); return a; } + +template N edge_of_triangle_with_angles(N alpha, N beta, N gamma) { + N of = (cos(alpha) + cos(beta) * cos(gamma)) / (sin(beta) * sin(gamma)); + if(hyperbolic) return acosh(of); + return acos(of); + } + +template N get_edgelen() { + N beta = get_deg(360)/valence; + return edge_of_triangle_with_angles (beta, get_deg(180)/edges, get_deg(180)/edges); + } + +template typename T::isometry cpush(int c, typename T::data::Number distance) { + return T::lorentz(c, T::data::Dim-1, distance); + } + +template struct cube_rotation_data_t { + std::vector> mapping; + }; + +template cube_rotation_data_t cube_rotation_data; + +template cube_rotation_data_t& build_cube_rotation() { + auto& crd = cube_rotation_data; + auto& crdm = crd.mapping; + // using N = typename T::data::Number; + if(crdm.empty()) { + crdm.emplace_back(hr::Id, T::id()); + for(int i=0; i(90)))); + if(is_new) crdm.emplace_back(U, T::apply(crdm[i].second, T::cspin90(j, k))); + } + if(isize(crdm) != 24) { + println(hlog, "the number of rotations found: ", isize(crdm)); + throw hr::hr_exception("wrong number of rotations"); + } + } + return crd; + } + +template U apply_move(cell *a, cell *b, U to) { + if(a == b) return to; + using N = typename T::data::Number; + if constexpr(test_dim == 4) { + auto& crdm = build_cube_rotation().mapping; + int ida = neighborId(a, b); + auto M = hr::currentmap->adj(a, ida); + for(int i0=0; i0(0, get_edgelen()), to); + to = T::apply(crdm[i0].second, to); + return to; + } + println(hlog, "tessf = ", hr::cgip->tessf); + println(hlog, "len = ", get_edgelen()); + throw hr::hr_exception("rotation not found"); + } + int ida = neighborId(a, b); + int idb = neighborId(b, a); + auto P1 = T::cspin(0, 1, idb * get_deg(360) / edges); + to = T::apply(P1, to); + auto P2 = cpush(0, get_edgelen()); + to = T::apply(P2, to); + auto P3 = T::cspin(1, 0, get_deg(180) + ida * get_deg(360) / edges); + to = T::apply(P3, to); + return to; + } + +template U apply_path(vector path, U to) { + for(int i=hr::isize(path)-2; i>=0; i--) + to = apply_move (path[i], path[i+1], to); + return to; + } + +template double test_sanity(int i) { + hr::indenter in(2); + + ld d1 = 1.25, d2 = 1.5; + + typename good::point gp = good::center(); + gp = good::apply(cpush(0, d1), gp); + gp = good::apply(cpush(1, d2), gp); + gp = good::apply(good::cspin(0, 1, get_deg(72)), gp); + + typename T::point p = T::center(); + p = T::apply(cpush(0, d1), p); + p = T::apply(cpush(1, d2), p); + p = T::apply(T::cspin(0, 1, get_deg(72)), p); + + double res = 0; + #define ADD(x, y) if(debug) println(hlog, "VS ", x, ",", y); res += pow( double(x-y), 2); + #define ADDA(x, y) if(debug) println(hlog, "VS ", x, ",", y); res += pow( cyclefix_on(double(x-y)), 2); + + if(debug) println(hlog, "p=", T::print(p)); + + ADD(T::get_coord(p, 0), good::get_coord(gp, 0)); + ADD(T::get_coord(p, 1), good::get_coord(gp, 1)); + ADD(T::get_coord(p, 2), good::get_coord(gp, 2)); + ADD(T::dist0(p), good::dist0(gp)); + ADDA(T::angle(p), good::angle(gp)); + + return res; + } + +template double test_consistency(int i) { + double res = 0; + using D = typename T::data; + + auto a = cpush(0, 1); + auto b = cpush(1, 1); + auto c = cpush(0, 1); + + auto s = T::apply(T::apply(a, b), c); + auto sp = T::apply(s, T::center()); + auto s1 = T::apply(a, T::apply(b, c)); + auto sp1 = T::apply(s1, T::center()); + auto sp2 = T::apply(a, T::apply(b, T::apply(c, T::center()))); + ADD(T::dist0(sp), T::dist0(sp1)); + ADD(T::dist0(sp), T::dist0(sp2)); + for(int i=0; i double test_tba(int id) { + std::mt19937 testr; testr.seed(id); + for(int i=0; itype; i++) { + vector p = {hr::cwt.at, hr::cwt.at->cmove(i), hr::cwt.at}; + auto h = apply_path(p, T::center()); + if(debug == 2) println(hlog, "i=", hr::lalign(3, i), " h = ", T::print(h), " distance =", T::dist0(h)); + if(T::dist0(h) >= 0 && T::dist0(h) < 0.1) continue; + exit(1); + } + return 999; + } + +template double test_loop_point(int id) { + std::mt19937 testr; testr.seed(id); + for(int i=0; i<100; i++) { + auto p1 = randomwalk(testr, hr::cwt.at, i); + auto p2 = random_return(testr, p1.back(), hr::cwt.at, 1/16., 1/32.); + auto p = vcon(p1, p2); + if(debug == 2) println(hlog, "path = ", p); + auto h = apply_path(vrev(p), T::center()); + if(debug == 2) println(hlog, "i=", hr::lalign(3, i), " h = ", T::print(h), " distance =", T::dist0(h)); + if(T::dist0(h) >= 0 && T::dist0(h) < 0.1) continue; + return i; + } + return 999; + } + +template double test_loop_iso(int id) { + std::mt19937 testr; testr.seed(id); + for(int i=0; i<100; i++) { + auto p1 = randomwalk(testr, hr::cwt.at, i); + auto p2 = random_return(testr, p1.back(), hr::cwt.at, 1/16., 1/32.); + auto p = vcon(p1, p2); + if(debug == 2) println(hlog, "path = ", p); + auto h = apply_path(vrev(p), T::id()); + auto hr = T::apply(h, T::center()); + // println(hlog, "i=", hr::lalign(3, i), " h=", hr); + if(debug == 2) println(hlog, "i=", hr::lalign(3, i), " hr = ", T::print(hr), " distance = ", T::dist0(hr)); + if(T::dist0(hr) >= 0 && T::dist0(hr) < 0.1) continue; + if(debug == 1) println(hlog, "i=", hr::lalign(3, i), " hr = ", T::print(hr), " distance = ", T::dist0(hr)); + return i; + } + return 999; + } + +template vector repeat_test(const F& f, int qty) { + vector res; + for(int i=0; i typename T::isometry random_rotation(std::mt19937& testr) { + using D = typename T::data; + using N = typename D::Number; + + if(D::Dim == 3) { + auto alpha = rand01(testr) * get_deg(360); + return T::cspin(0, 1, alpha); + } + + auto x = T::id(); + for(int i=0; i<100; i++) { + int c0 = testr() % (D::Dim-1); + int c1 = testr() % (D::Dim-1); + if(c0 == c1) continue; + auto len = rand01(testr) * get_deg(360); + x = T::apply(T::cspin(c0, c1, len), x); + } + return x; + } + +template typename T::isometry random_iso(std::mt19937& testr) { + auto x = T::id(); + using D = typename T::data; + using N = typename D::Number; + for(int i=0; i<100; i++) { + int c0 = testr() % D::Dim; + int c1 = testr() % D::Dim; + if(c0 == c1) continue; + if(c0 == D::Flipped) std::swap(c0, c1); + N len = c1 < D::Flipped ? rand01(testr) * get_deg(360) : (rand01(testr)-N(0.5)) * N(0.25); + if(c1 == D::Flipped) x = T::apply(T::lorentz(c0, c1, len), x); + else x = T::apply(T::cspin(c0, c1, len), x); + } + return x; + } + + +template std::string test_count(int id) { + std::mt19937 testr; testr.seed(id); + hr::shstream out; + auto A = random_iso(testr); + auto B = random_iso(testr); + auto C = random_iso(testr); + auto P = T::apply(C, T::center()); + for(int i=0; i<9; i++) { + counts.clear(); + for(auto& i: cbc) i = 0; + std::string s; + switch(i) { + case 0: s = "spin"; T::cspin(0, 1, countfloat(.5)); break; + case 1: s = "L0"; T::lorentz(0, T::data::Dim-1, countfloat(.5)); break; + case 2: s = "L1"; T::lorentz(1, T::data::Dim-1, countfloat(.5)); break; + case 3: s = "ip"; T::apply(A, P); break; + case 4: s = "ii"; T::apply(A, B); break; + case 5: s = "d0"; T::dist0(P); break; + case 6: s = "angle"; T::angle(P); break; + case 7: s = "inverse"; T::inverse(A); break; + case 8: s = "push"; T::push(P); break; + } + if(i) print(out, " "); + if(1) { + print(out, s, "("); + bool nsp = false; + for(int i=1; i<5; i++) if(cbc[i]) { + if(nsp) print(out, " "); + print(out, cbc[i], hr::s0+".AMDF"[i]); + nsp = true; + } + print(out, ")"); + } + } + return out.s; + } + +template bool closeto(A a, B b) { return abs(a-b) < 0.1; } + +template bool closeto_angle(A a, B b) { return abs(cyclefix_on(double(a-b))) < 0.1; } + +template double test_angledist(int id) { + std::mt19937 testr; testr.seed(id); + for(int i=1; i<1000; i += (1 + i/5)) { + + auto p = randomwalk(testr, hr::cwt.at, i); + auto h = apply_path(vrev(p), T::center()); + + auto gh = apply_path(vrev(p), good::center()); + + if(debug == 2) { + println(hlog, "good: ", good::print(gh), " dist = ", good::dist0(gh), " angle = ", good::angle(gh)); + println(hlog, "test: ", T::print(h), " dist = ", T::dist0(h), " angle = ", T::angle(h), " [i=", i, "]"); + } + + if(closeto(good::dist0(gh), T::dist0(h)) && closeto_angle(good::angle(gh), T::angle(h))) continue; + + if(debug == 1) { + println(hlog, "good: ", good::print(gh), " dist = ", good::dist0(gh), " angle = ", good::angle(gh)); + println(hlog, "test: ", T::print(h), " dist = ", T::dist0(h), " angle = ", T::angle(h), " [i=", i, "]"); + } + + return i; + } + return 999; + } + +#define TEST_VARIANTS(x,D,q,t,rn, r) \ + nm = nmInvariant; println(hlog, rn, "invariant: ", repeat_test(x>, q)); \ + nm = nmForced; println(hlog, rn, "forced : ", repeat_test(x>, q)); \ + nm = nmWeak; println(hlog, rn, "weak : ", repeat_test(x>, q)); \ + nm = nmFlatten; println(hlog, rn, "flatten : ", repeat_test(x>, q)); \ + nm = nmCareless; println(hlog, rn, "careless : ", repeat_test(x>, q)); \ + nm = nmBinary; println(hlog, rn, "binary : ", repeat_test(x>, q)); + +/* +#define TEST_ALL(x,D,q,t) \ + println(hlog, "HyperRogue: ", repeat_test(x>, q)); \ + polar_mod = polar_choose = false; println(hlog, "high polar: ", repeat_test(x>, q)); \ + if(test_dim == 3) { polar_mod = polar_choose = false; println(hlog, "low polar : ", repeat_test(x>, q)); } +*/ + +// println(hlog, "HyperRogue : ", repeat_test(x>, q)); + +#define TEST_ALL(x,D,q,t) \ + fix_matrices = true; TEST_VARIANTS(x,D,q,t,"linear+F ", rep_linear) \ + fix_matrices = false; TEST_VARIANTS(x,D,q,t,"linear-F ", rep_linear) \ + TEST_VARIANTS(x,D,q,t,"mixed ", rep_mixed) \ + TEST_VARIANTS(x,D,q,t,"Clifford ", rep_clifford) \ + nm = nmFlatten; println(hlog, "Clifford gyro : ", repeat_test(x>, q)); \ + nm = nmInvariant; println(hlog, "halfplane invariant: ", repeat_test(x>, q)); \ + polar_choose = false; println(hlog, "polar basic : ", repeat_test(x>, q)); \ + polar_choose = true; println(hlog, "polar improved : ", repeat_test(x>, q)); \ + if(test_dim == 3) { \ + polar_mod = false; polar_choose = false; println(hlog, "polar F/F : ", repeat_test(x>, q)); \ + polar_mod = false; polar_choose = true; println(hlog, "polar F/T : ", repeat_test(x>, q)); \ + polar_mod = true; polar_choose = false; println(hlog, "polar T/F : ", repeat_test(x>, q)); \ + polar_mod = true; polar_choose = true; println(hlog, "polar T/T : ", repeat_test(x>, q)); \ + } + +template double test_distances(int id, int a) { + std::mt19937 testr; testr.seed(id); + using N = typename T::data::Number; + + for(int i=1; i<1000; i ++) { + + auto R = random_rotation(testr); + auto dif = exp(N(-1) * i) + get_deg(a); + + auto p1 = T::apply(T::apply(R, cpush(0, N(i))), T::center()); + auto p2 = T::apply(T::apply(R, T::apply(T::cspin(0, 1, dif), cpush(0, N(i)))), T::center()); + auto pd = T::apply(T::inverse(T::push(p1)), p2); + auto d = T::dist0(pd); + + // for good we do not need R actually + auto gp1 = good::apply(cpush(0, N(i)), good::center()); + auto gp2 = good::apply(good::apply(good::cspin(0, 1, dif), cpush(0, N(i))), good::center()); + auto gd = good::dist0(good::apply(good::inverse(good::push(gp1)), gp2)); + + if(debug == 2) println(hlog, T::print(p1), " ... ", T::print(p2), " = ", T::print(pd), " d=", d, " [i=", i, " dif=", dif, "]"); + + if(closeto(d, gd)) continue; + + return i; + } + return 999; + } + +template double test_similarity(int id) { return test_distances(id, 0); } +template double test_dissimilarity(int id) { return test_distances(id, 180); } +template double test_other(int id) { return test_distances(id, 1); } + +template double test_walk(int id) { + std::mt19937 testr; testr.seed(id); + + ld step = 1/16.; + // mover-relative to cell-relative + auto R0 = random_rotation(testr); + cell *c0 = hr::cwt.at; + auto R1 = T::apply(R0, cpush(0, step/2)); + cell *c1 = hr::cwt.at; + + int i = 0; + int lastchange = 0; + while(i< lastchange + 1000 && i < 10000 && celldistance(c0, c1) < 3) { + // println(hlog, "iteration ", i, " in ", c0, " vs ", c1); + auto rebase = [&] (typename T::isometry& R, cell*& c, int id) { + ld d = T::dist0(T::apply(R, T::center())); + for(int dir=0; dirtype; dir++) { + cell *altc = c->cmove(dir); + auto altR = apply_move(altc, c, R); + ld altd = T::dist0(T::apply(altR, T::center())); + if(altd < d + 1/256.) { + R = altR; c = altc; lastchange = i; return; + } + } + }; + R0 = T::apply(R0, cpush(0, step)); rebase(R0, c0, 0); + R1 = T::apply(R1, cpush(0, step)); rebase(R1, c1, 1); + i++; + } + + return i; + } + +template double test_close(int id) { + std::mt19937 testr; testr.seed(id); + + cell *c = hr::cwt.at; + int phase = 0; + auto p0 = T::apply(cpush(0, 1/8.), T::center()); + auto p = p0; + + int steps = 0; + const int maxdist = id + 1; + int errors = 0; + + while(steps < 10000) { + int d = testr() % c->type; + cell *c1 = c->cmove(d); + + bool do_move = false; + + switch(phase) { + case 0: + /* always move */ + do_move = true; + if(celldistance(c1, hr::cwt.at) == maxdist) phase = 1; + break; + + case 1: + /* move only towards the center */ + int d0 = celldistance(c, hr::cwt.at); + int d1 = celldistance(c1, hr::cwt.at); + do_move = d1 < d0; + if(d1 == 0) phase = 0; + break; + } + + if(do_move) { + p = apply_move(c1, c, p); c = c1; steps++; + if(debug == 2) println(hlog, "dist = ", celldistance(c, hr::cwt.at), " dist = ", T::dist0(p)); + if(c == hr::cwt.at) { + auto d = T::dist0(p); + auto a = T::angle(p); + if(!closeto(d, 1/8.) || !closeto_angle(a, 0)) { + errors++; phase = 0; c = hr::cwt.at; p = p0; + } + } + } + } + + return errors; + } + +void run_all_tests() { + prepare_tests(); + + // println(hlog, "test_sanity"); TEST_ALL(test_sanity, data, 1, ld); + + // println(hlog, "test_consistency"); TEST_ALL(test_consistency, data, 1, ld); + + println(hlog, "test_loop_iso"); TEST_ALL(test_loop_iso, data, 20, int); + + println(hlog, "test_loop_point"); TEST_ALL(test_loop_point, data, 20, int); + + println(hlog, "test_angledist"); TEST_ALL(test_angledist, data, 3, int); + + println(hlog, "test_similarity"); TEST_ALL(test_similarity, data, 20, int); + + println(hlog, "test_dissimilarity"); TEST_ALL(test_dissimilarity, data, 20, int); + + println(hlog, "test_other"); TEST_ALL(test_other, data, 20, int); + + println(hlog, "test_walk"); TEST_ALL(test_walk, data, 20, int); + + println(hlog, "test_close"); TEST_ALL(test_close, data, 20, int); + + println(hlog, "test_count"); TEST_ALL(test_count, countdata, 1, std::string); + } + +}