Commits

Barry Schwartz committed fc4d0eb

A changed implementation of Pure multivector addition. I haven’t bothered to do timings, but it is meant more as a prototype of how to do these things.

Comments (0)

Files changed (5)

 	_dst->c[1] = a->c[1];
 
 }
+void add_bvec2_bvec2(bvec2 *_dst, const bvec2 *a, const bvec2 *b)
+{
+	_dst->c[0] = (a->c[0]+b->c[0]);
+
+}
 void subtract_mvec2_mvec2(mvec2 *_dst, const mvec2 *a, const mvec2 *b)
 {
 	_dst->c[0] = (a->c[0]-b->c[0]);
 	_dst->c[1] = a->c[1];
 
 }
+void subtract_bvec2_bvec2(bvec2 *_dst, const bvec2 *a, const bvec2 *b)
+{
+	_dst->c[0] = (a->c[0]-b->c[0]);
+
+}
 void gp_mvec3c_mvec3c(mvec3c *_dst, const mvec3c *a, const mvec3c *b)
 {
 	double c[32];
  */
 void add_rotor2_double(rotor2 *_dst, const rotor2 *a, const double b);
 /**
+ * Returns bvec2 + bvec2.
+ */
+void add_bvec2_bvec2(bvec2 *_dst, const bvec2 *a, const bvec2 *b);
+/**
  * Returns mvec2 - mvec2.
  */
 void subtract_mvec2_mvec2(mvec2 *_dst, const mvec2 *a, const mvec2 *b);
  */
 void subtract_rotor2_double(rotor2 *_dst, const rotor2 *a, const double b);
 /**
+ * Returns bvec2 - bvec2.
+ */
+void subtract_bvec2_bvec2(bvec2 *_dst, const bvec2 *a, const bvec2 *b);
+/**
  * Returns geometric product of mvec3c and mvec3c.
  */
 void gp_mvec3c_mvec3c(mvec3c *_dst, const mvec3c *a, const mvec3c *b);
   <function name="add" arg1="mvec2" arg2="double"/>
   <function name="add" arg1="double" arg2="rotor2"/>
   <function name="add" arg1="rotor2" arg2="double"/>
+  <function name="add" arg1="bvec2" arg2="bvec2"/>
 
   <function name="subtract" arg1="mvec2" arg2="mvec2"/>
   <function name="subtract" arg1="vec2" arg2="vec2"/>
   <function name="subtract" arg1="mvec2" arg2="double"/>
   <function name="subtract" arg1="double" arg2="rotor2"/>
   <function name="subtract" arg1="rotor2" arg2="double"/>
+  <function name="subtract" arg1="bvec2" arg2="bvec2"/>
 
   <function name="gp" arg1="mvec3c" arg2="mvec3c"/>
   <function name="gp" arg1="mvec2" arg2="mvec2"/>

pure/ebeno/geometric_algebra.pure

 using ebeno::geometric_algebra_h;
 using ebeno::geometric_algebra_1;
 using math;
-using system;
 
 namespace ebeno;
 

pure/generate_ga_code.pure

 
 //-------------------------------------------------------------------------
 
-// FIXME: Consider calling the C code instead of compiling ‘double’
-// versions.
+// A function that creates multivector addition functions.
 
-// A function that creates multivector addition functions.
 make_multivector_addition type1 type2 indices1::list indices2::list type_out =
-  sprintf "%s {%s} ::+ %s {%s} |\n"
-  (str type1, a_elemlist_double, str type2, b_elemlist_double) +
-  sprintf "%s {%s} ::+ %s {%s} = %s {%s};\n"
-  (str type1, a_elemlist, str type2, b_elemlist, str type_out, elemlist_out)
+  (if type1 === type2 then
+     str type1 + " c1@{" + a_elemlist + "} ::+ " +
+     str type2 + " c2@{" + b_elemlist + "} =\n"
+   else
+     str type1 + " c1@{" + a_elemlist + "} ::+ " +
+     str type2 + " c2@{" + b_elemlist + "} |\n" +
+     str type2 + " c2@{" + b_elemlist + "} ::+ " +
+     str type1 + " c1@{" + a_elemlist + "} =\n") +
+  "  if dmatrixp c1 && dmatrixp c2 then\n" +
+  "    (" + ebeno_func + " (pointer c, pointer cp1, pointer cp2) $$ " + str type_out + " c\n" +
+  "     when\n" +
+  "       c::dmatrix = rowvector $ repeatn " + str (m + 1) + " 0.0;\n" +
+  "       cp1::dmatrix = pack c1;\n" +
+  "       cp2::dmatrix = pack c2;\n" +
+  "     end)\n" +
+  "  else\n" +
+  "    " + str type_out + " {" + elemlist_out + "}\n" +
+  "  if ~cmatrixp c1 && ~cmatrixp c2;\n"
   when
+    ebeno_func = "add_" + str type1 + "_" + str type2;
     a_elem = ["a" + str i | i = indices1];
     b_elem = ["b" + str i | i = indices2];
     a_elemlist = foldl (\el s -> el + ("," + s)) (a_elem!0) (drop 1 a_elem);
       end;
     elemlist_out = foldl (\el s -> el + ("," + s)) (elem_out!0) (drop 1 elem_out);
   end;
+
 make_multivector_addition type1 type2 indices1::list indices2 type_out =
-  sprintf "%s {%s} ::+ r::double |\n" (str type1, a_elemlist_double) +
-  sprintf "%s {%s} ::+ r::real = %s {%s};\n" (str type1, a_elemlist, str type_out, elemlist_out)
+  str type1 + " c1@{" + a_elemlist + "} ::+ r::real |\n" +
+  "r::real ::+ " + str type1 + " c1@{" + a_elemlist + "} =\n" +
+  "  if dmatrixp c1 then\n" +
+  "    (" + ebeno_func + " (pointer c, pointer cp1, double r) $$ " + str type_out + " c\n" +
+  "     when\n" +
+  "       c::dmatrix = rowvector $ repeatn (#c1) 0.0;\n" +
+  "       cp1::dmatrix = pack c1;\n" +
+  "     end)\n" +
+  "  else\n" +
+  "    " + str type_out + " {" + elemlist_out + "}\n" +
+  "  if ~cmatrixp c1;\n"
   when
+    ebeno_func = "add_" + str type1 + "_double";
     a_elem = ["a" + str i | i = indices1];
     a_elemlist = foldl (\el s -> el + ("," + s)) (a_elem!0) (drop 1 a_elem);
     a_elemlist_double = foldl (\el s -> el + ("," + s + "::double")) (a_elem!0 + "::double") (drop 1 a_elem);
     elemlist_out = foldl (\el s -> el + ("," + s)) (a_elem!0 + " + r") (drop 1 a_elem);
   end;
-make_multivector_addition type1 type2 indices1 indices2::list type_out =
-  sprintf "r::double ::+ %s {%s} |\n" (str type2, b_elemlist_double) +
-  sprintf "r::real ::+ %s {%s} = %s {%s};\n" (str type2, b_elemlist, str type_out, elemlist_out)
-  when
-    b_elem = ["b" + str i | i = indices2];
-    b_elemlist = foldl (\el s -> el + ("," + s)) (b_elem!0) (drop 1 b_elem);
-    b_elemlist_double = foldl (\el s -> el + ("," + s + "::double")) (b_elem!0 + "::double") (drop 1 b_elem);
-    elemlist_out = foldl (\el s -> el + ("," + s)) ("r + " + b_elem!0) (drop 1 b_elem);
-  end;
 
 //-------------------------------------------------------------------------
 
 // versions.
 
 // A function that creates multivector subtraction functions.
+
 make_multivector_subtraction type1 type2 indices1::list indices2::list type_out =
   sprintf "%s {%s} ::- %s {%s} |\n"
   (str type1, a_elemlist_double, str type2, b_elemlist_double) +
 puts $ make_multivector_addition mvec2 vec2 [0,1,2,3] [1,2] mvec2;
 puts $ make_multivector_addition mvec2 bvec2 [0,1,2,3] [3] mvec2;
 puts $ make_multivector_addition mvec2 rotor2 [0,1,2,3] [0,3] mvec2;
-
-puts $ make_multivector_addition vec2 mvec2 [1,2] [0,1,2,3] mvec2;
 puts $ make_multivector_addition vec2 vec2 [0,1] [0,1] vec2;
-
-puts $ make_multivector_addition bvec2 mvec2 [3] [0,1,2,3] mvec2;
 puts $ make_multivector_addition bvec2 bvec2 [0] [0] bvec2;
 puts $ make_multivector_addition bvec2 rotor2 [1] [0,1] rotor2;
-
-puts $ make_multivector_addition rotor2 mvec2 [0,3] [0,1,2,3] mvec2;
-puts $ make_multivector_addition rotor2 bvec2 [0,1] [1] rotor2;
 puts $ make_multivector_addition rotor2 rotor2 [0,1] [0,1] rotor2;
 
 puts $ make_multivector_addition mvec2 real [0,1,2,3] () mvec2;