From 469fdea0a118e64b941840321e02fdf56e79263f Mon Sep 17 00:00:00 2001
From: Bill Currie <bill@taniwha.org>
Date: Wed, 23 Aug 2023 21:52:33 +0900
Subject: [PATCH] [qfcc] Implement 2d PGA outer products

The specialization is really just in the layout.
---
 tools/qfcc/source/expr_algebra.c | 39 ++++++++++++++++++++++++++++++++
 tools/qfcc/test/algtypes.r       |  9 ++++++++
 2 files changed, 48 insertions(+)

diff --git a/tools/qfcc/source/expr_algebra.c b/tools/qfcc/source/expr_algebra.c
index 04aa61fd6..de8f5d876 100644
--- a/tools/qfcc/source/expr_algebra.c
+++ b/tools/qfcc/source/expr_algebra.c
@@ -422,6 +422,40 @@ static void (*pga3_wedge_funcs[6][6])(expr_t**,expr_t*,expr_t*,algebra_t*) = {
 	},
 };
 
+// vector-bivector wedge is commutative
+#define pga2_x_y_w_wedge_yw_wx_xy pga2_yw_wx_xy_wedge_x_y_w
+static void
+pga2_yw_wx_xy_wedge_x_y_w (expr_t **c, expr_t *a, expr_t *b, algebra_t *alg)
+{
+	c[3] = dot_expr (algebra_mvec_type (alg, 0x08), a, b);
+}
+
+static void
+pga2_x_y_w_wedge_x_y_w (expr_t **c, expr_t *a, expr_t *b, algebra_t *alg)
+{
+	auto wedge_type = algebra_mvec_type (alg, 0x01);
+	c[0] = new_binary_expr (CROSS, a, b);
+	c[0]->e.expr.type = wedge_type;
+}
+
+static void (*pga2_wedge_funcs[4][4])(expr_t**,expr_t*,expr_t*,algebra_t*) = {
+	[0] = {
+		[1] = scale_component,
+		[2] = pga2_yw_wx_xy_wedge_x_y_w,
+	},
+	[1] = {
+		[0] = scale_component,
+		[1] = scale_component,
+		[2] = scale_component,
+		[3] = scale_component,
+	},
+	[2] = {
+		[0] = pga2_x_y_w_wedge_yw_wx_xy,
+		[1] = scale_component,
+		[2] = pga2_x_y_w_wedge_x_y_w,
+	},
+};
+
 static void
 component_wedge (expr_t **c, expr_t *a, expr_t *b, algebra_t *algebra)
 {
@@ -436,6 +470,11 @@ component_wedge (expr_t **c, expr_t *a, expr_t *b, algebra_t *algebra)
 			pga3_wedge_funcs[ga][gb] (c, a, b, algebra);
 		}
 	} else if (p == 2 && m == 0 && z == 1) {
+		int ga = get_group (get_type (a), algebra);
+		int gb = get_group (get_type (b), algebra);
+		if (pga2_wedge_funcs[ga][gb]) {
+			pga2_wedge_funcs[ga][gb] (c, a, b, algebra);
+		}
 	} else {
 	}
 }
diff --git a/tools/qfcc/test/algtypes.r b/tools/qfcc/test/algtypes.r
index a03395a22..3e32c7b7b 100644
--- a/tools/qfcc/test/algtypes.r
+++ b/tools/qfcc/test/algtypes.r
@@ -9,6 +9,9 @@
 typedef @algebra(float(3,0,1)) PGA;
 typedef @algebra(float(4,1)) CGA;
 
+typedef @algebra(double(2,0,1)) PGA2;
+PGA2 pga2;
+
 float sin(float x) = #0;
 
 int
@@ -34,5 +37,11 @@ main (void)
 		auto mvec = rx + ry;
 #endif
 	}
+	@algebra (PGA2) {
+		auto l1 = e1 + 2 * e2 + 5 * e0;
+		auto l2 = 3 * e1 - e2 + 10 * e0;
+		auto p = l1∧l2;
+		pga2 = p + (1 + p)∧l1;
+	}
 	return 0;		// to survive and prevail :)
 }