From 40f04ff33e4ebe0e555933096a1f6f0a0b5da6ca Mon Sep 17 00:00:00 2001 From: Bill Currie Date: Fri, 9 Feb 2024 13:49:36 +0900 Subject: [PATCH] [qfcc] Implement undual and use for regressive product I'm not sure the regressive product is right (overall sign), but that's actually partly a problem in the math itself (duals and the regressive product still get poked at, so it may be just a matter of interpretation). --- tools/qfcc/include/algebra.h | 1 + tools/qfcc/source/expr.c | 2 + tools/qfcc/source/expr_algebra.c | 38 +++++++++--- tools/qfcc/source/qc-lex.l | 2 + tools/qfcc/source/qc-parse.y | 3 +- tools/qfcc/test/pga2d.r | 102 ++++++++++++++++++++++++++++++- tools/qfcc/test/pga3d.r | 75 +++++++++++++++++++++++ 7 files changed, 212 insertions(+), 11 deletions(-) diff --git a/tools/qfcc/include/algebra.h b/tools/qfcc/include/algebra.h index 152f58234..032c87cf5 100644 --- a/tools/qfcc/include/algebra.h +++ b/tools/qfcc/include/algebra.h @@ -115,6 +115,7 @@ pr_uint_t get_group_mask (const type_t *type, algebra_t *algebra) __attribute__( const expr_t *algebra_binary_expr (int op, const expr_t *e1, const expr_t *e2); const expr_t *algebra_negate (const expr_t *e); const expr_t *algebra_dual (const expr_t *e); +const expr_t *algebra_undual (const expr_t *e); const expr_t *algebra_reverse (const expr_t *e); const expr_t *algebra_cast_expr (type_t *dstType, const expr_t *e); const expr_t *algebra_assign_expr (const expr_t *dst, const expr_t *src); diff --git a/tools/qfcc/source/expr.c b/tools/qfcc/source/expr.c index 5ae429c9e..663c4c228 100644 --- a/tools/qfcc/source/expr.c +++ b/tools/qfcc/source/expr.c @@ -2058,6 +2058,8 @@ bitnot_expr: return algebra_reverse (e); case DUAL: return algebra_dual (e); + case UNDUAL: + return algebra_undual (e); } internal_error (e, 0); } diff --git a/tools/qfcc/source/expr_algebra.c b/tools/qfcc/source/expr_algebra.c index 853a2b2c9..67e8fa9cb 100644 --- a/tools/qfcc/source/expr_algebra.c +++ b/tools/qfcc/source/expr_algebra.c @@ -2596,7 +2596,7 @@ regressive_product (const expr_t *e1, const expr_t *e2) auto a = algebra_dual (e1); auto b = algebra_dual (e2); auto c = outer_product (a, b); - return algebra_dual (c); + return algebra_undual (c); } static const expr_t * @@ -2762,12 +2762,13 @@ algebra_negate (const expr_t *e) return mvec_gather (n, algebra); } -const expr_t * -algebra_dual (const expr_t *e) +static const expr_t * +hodge_dual (const expr_t *e, bool undual) { auto algebra = algebra_context (get_type (e)); if (!algebra) { - return error (e, "cannot take the dual of a scalar without context"); + return error (e, "cannot take the %s of a scalar without context", + undual ? "undual" : "dual"); } auto layout = &algebra->layout; @@ -2776,6 +2777,7 @@ algebra_dual (const expr_t *e) e = mvec_expr (e, algebra); mvec_scatter (a, e, algebra); + int dim = algebra->dimension - algebra->zero; pr_uint_t I_mask = (1u << algebra->dimension) - 1; for (int i = 0; i < layout->count; i++) { if (!a[i]) { @@ -2784,21 +2786,39 @@ algebra_dual (const expr_t *e) auto group = &layout->groups[i]; //FIXME assumes groups are mono-grade (either come up with something //or reject mixed-grade groups) - pr_uint_t mask = I_mask ^ group->blades[0].mask; - int dual_ind = layout->group_map[layout->mask_map[mask]][0]; + auto blade = group->blades[0]; + pr_uint_t d_mask = I_mask ^ blade.mask; + int dual_ind = layout->group_map[layout->mask_map[d_mask]][0]; auto dual_group = &layout->groups[dual_ind]; auto dual_type = algebra_mvec_type (algebra, dual_group->group_mask); auto dual = cast_expr (dual_type, a[i]); - if (algebra_count_flips (algebra, group->blades[0].mask, mask) & 1) { - b[dual_ind] = neg_expr (dual); + pr_uint_t mask = algebra_count_flips (algebra, blade.mask, d_mask); + if (mask & 1) { + dual = neg_expr (dual); } else { - b[dual_ind] = dual; + dual = dual; } + if (undual && ((dim * algebra_get_grade (get_type (a[i]))) & 1)) { + dual = neg_expr (dual); + } + b[dual_ind] = dual; } return mvec_gather (b, algebra); } +const expr_t * +algebra_dual (const expr_t *e) +{ + return hodge_dual (e, false); +} + +const expr_t * +algebra_undual (const expr_t *e) +{ + return hodge_dual (e, true); +} + static void set_sign (pr_type_t *val, int sign, const type_t *type) { diff --git a/tools/qfcc/source/qc-lex.l b/tools/qfcc/source/qc-lex.l index e99e3a557..f312a2146 100644 --- a/tools/qfcc/source/qc-lex.l +++ b/tools/qfcc/source/qc-lex.l @@ -443,6 +443,8 @@ static keyword_t qf_keywords[] = { {"@regressive", REGRESSIVE, }, {"@geometric", GEOMETRIC, }, {"@algebra", ALGEBRA, }, + {"@dual", DUAL, }, + {"@undual", UNDUAL, }, }; // These keywors are always available. Other than the @ keywords, they diff --git a/tools/qfcc/source/qc-parse.y b/tools/qfcc/source/qc-parse.y index aa4011aff..89af3eb0b 100644 --- a/tools/qfcc/source/qc-parse.y +++ b/tools/qfcc/source/qc-parse.y @@ -145,7 +145,7 @@ int yylex (void); %left '+' '-' %left '*' '/' '%' MOD SCALE GEOMETRIC QMUL QVMUL VQMUL %left HADAMARD CROSS DOT WEDGE REGRESSIVE -%right SIZEOF UNARY INCOP REVERSE STAR DUAL +%right SIZEOF UNARY INCOP REVERSE STAR DUAL UNDUAL %left HYPERUNARY %left '.' '(' '[' @@ -1730,6 +1730,7 @@ unary_expr | unary_expr INCOP { $$ = incop_expr ($2, $1, 1); } | unary_expr REVERSE { $$ = unary_expr (REVERSE, $1); } | DUAL cast_expr %prec UNARY { $$ = unary_expr (DUAL, $2); } + | UNDUAL cast_expr %prec UNARY { $$ = unary_expr (UNDUAL, $2); } | '+' cast_expr %prec UNARY { $$ = $2; } | '-' cast_expr %prec UNARY { $$ = unary_expr ('-', $2); } | '!' cast_expr %prec UNARY { $$ = unary_expr ('!', $2); } diff --git a/tools/qfcc/test/pga2d.r b/tools/qfcc/test/pga2d.r index 09f869090..527b08f1d 100644 --- a/tools/qfcc/test/pga2d.r +++ b/tools/qfcc/test/pga2d.r @@ -20,6 +20,96 @@ typedef union { }; } oddgrades_t; +static int +test_dual (void) +{ + @algebra (PGA) { + auto a = 1f; + auto a0 = e0; + auto a1 = e1; + auto a2 = e2; + auto a12 = e12; + auto a01 = e01; + auto a20 = e20; + auto a012 = e012; + #define TEST_DUAL(x, y) \ + if (⋆x != y) { \ + printf ("⋆" #x " != " #y "\n"); \ + return 1; \ + } + TEST_DUAL (a, e012); + TEST_DUAL (a0, e12); + TEST_DUAL (a1, e20); + TEST_DUAL (a2, e01); + TEST_DUAL (a20, e1); + TEST_DUAL (a01, e2); + TEST_DUAL (a12, e0); + TEST_DUAL (a012, 1); + #undef TEST_DUAL + + #define TEST_DUAL(x) \ + if (x * ⋆x != e012) { \ + printf (#x " * ⋆" #x " != e012\n"); \ + return 1; \ + } + TEST_DUAL (a); + TEST_DUAL (a0); + TEST_DUAL (a1); + TEST_DUAL (a2); + TEST_DUAL (a20); + TEST_DUAL (a01); + TEST_DUAL (a12); + TEST_DUAL (a012); + #undef TEST_DUAL + } + return 0; +} + +static int +test_undual (void) +{ + @algebra (PGA) { + auto a = 1f; + auto a0 = e0; + auto a1 = e1; + auto a2 = e2; + auto a12 = e12; + auto a01 = e01; + auto a20 = e20; + auto a012 = e012; + #define TEST_UNDUAL(x, y) \ + if (@undual x != y) { \ + printf ("@undual " #x " != " #y "\n"); \ + return 1; \ + } + TEST_UNDUAL (a, e012); + TEST_UNDUAL (a0, e12); + TEST_UNDUAL (a1, e20); + TEST_UNDUAL (a2, e01); + TEST_UNDUAL (a20, e1); + TEST_UNDUAL (a01, e2); + TEST_UNDUAL (a12, e0); + TEST_UNDUAL (a012, 1); + #undef TEST_UNDUAL + + #define TEST_UNDUAL(x) \ + if (@undual x * x != e012) { \ + printf ("@undual " #x " * " #x " != e012\n"); \ + return 1; \ + } + TEST_UNDUAL (a); + TEST_UNDUAL (a0); + TEST_UNDUAL (a1); + TEST_UNDUAL (a2); + TEST_UNDUAL (a20); + TEST_UNDUAL (a01); + TEST_UNDUAL (a12); + TEST_UNDUAL (a012); + #undef TEST_DUAL + } + return 0; +} + int main (void) { @@ -269,7 +359,17 @@ main (void) } auto line = bvec ∨ bvecb; if ((dvec3)line != '-11 8 34'd) { - printf ("bvec ∨ bvecb != '-11 8 34': %lv\n", line); + printf ("bvec(%lv) ∨ bvecb(%lv) != '-11 8 34': %lv\n", bvec, bvecb, line); + return 1; + } + + if (test_dual ()) { + printf ("dual failed\n"); + return 1; + } + + if (test_undual ()) { + printf ("dual failed\n"); return 1; } return 0; // to survive and prevail :) diff --git a/tools/qfcc/test/pga3d.r b/tools/qfcc/test/pga3d.r index 42031986c..81c593173 100644 --- a/tools/qfcc/test/pga3d.r +++ b/tools/qfcc/test/pga3d.r @@ -434,6 +434,76 @@ test_dual (void) TEST_DUAL (a021); TEST_DUAL (a123); TEST_DUAL (a0123); + #undef TEST_DUAL + } + return 0; +} + +static int +test_undual (void) +{ + @algebra (PGA) { + auto a = 1f; + auto a0 = e0; + auto a1 = e1; + auto a2 = e2; + auto a3 = e3; + auto a23 = e23; + auto a31 = e31; + auto a12 = e12; + auto a01 = e01; + auto a02 = e02; + auto a03 = e03; + auto a032 = e032; + auto a013 = e013; + auto a021 = e021; + auto a123 = e123; + auto a0123 = e0123; + #define TEST_UNDUAL(x, y) \ + if (@undual x != y) { \ + printf ("@undual" #x " != " #y "\n"); \ + return 1; \ + } + TEST_UNDUAL (a, e0123); + TEST_UNDUAL (a0, -e123); + TEST_UNDUAL (a1, -e032); + TEST_UNDUAL (a2, -e013); + TEST_UNDUAL (a3, -e021); + TEST_UNDUAL (a23, e01); + TEST_UNDUAL (a31, e02); + TEST_UNDUAL (a12, e03); + TEST_UNDUAL (a01, e23); + TEST_UNDUAL (a02, e31); + TEST_UNDUAL (a03, e12); + TEST_UNDUAL (a032, e1); + TEST_UNDUAL (a013, e2); + TEST_UNDUAL (a021, e3); + TEST_UNDUAL (a123, e0); + TEST_UNDUAL (a0123, 1); + #undef TEST_UNDUAL + + #define TEST_UNDUAL(x) \ + if (@undual x * x != e0123) { \ + printf ("@undual " #x " * " #x " != e0123\n"); \ + return 1; \ + } + TEST_UNDUAL (a); + TEST_UNDUAL (a0); + TEST_UNDUAL (a1); + TEST_UNDUAL (a2); + TEST_UNDUAL (a3); + TEST_UNDUAL (a23); + TEST_UNDUAL (a31); + TEST_UNDUAL (a12); + TEST_UNDUAL (a01); + TEST_UNDUAL (a02); + TEST_UNDUAL (a03); + TEST_UNDUAL (a032); + TEST_UNDUAL (a013); + TEST_UNDUAL (a021); + TEST_UNDUAL (a123); + TEST_UNDUAL (a0123); + #undef TEST_UNDUAL } return 0; } @@ -549,6 +619,11 @@ main (void) return 1; } + if (test_undual ()) { + printf ("dual failed\n"); + return 1; + } + if (test_basics ()) { printf ("basics failed\n"); return 1;