Index: libcfa/src/vec/vec.hfa
===================================================================
--- libcfa/src/vec/vec.hfa	(revision 3376ec9679544e4a3573e3f6ba122456afec811f)
+++ libcfa/src/vec/vec.hfa	(revision 3376ec9679544e4a3573e3f6ba122456afec811f)
@@ -0,0 +1,120 @@
+#pragma once
+
+#include <math.hfa>
+
+trait fromint(otype T) {
+    void ?{}(T&, int);
+};
+trait zeroinit(otype T) {
+    void ?{}(T&, zero_t);
+};
+trait zero_assign(otype T) {
+    T ?=?(T&, zero_t);
+};
+trait subtract(otype T) {
+    T ?-?(T, T);
+};
+trait negate(otype T) {
+    T -?(T);
+};
+trait add(otype T) {
+    T ?+?(T, T);
+};
+trait multiply(otype T) {
+    T ?*?(T, T);
+};
+trait divide(otype T) {
+    T ?/?(T, T);
+};
+trait lessthan(otype T) {
+    int ?<?(T, T);
+};
+trait equality(otype T) {
+    int ?==?(T, T);
+};
+trait sqrt(otype T) {
+    T sqrt(T);
+};
+
+static inline {
+// int
+int ?=?(int& n, zero_t) { return n = 0.f; }
+/* float */
+void ?{}(float& a, int b) { a = b; }
+float ?=?(float& n, zero_t) { return n = 0.f; }
+/* double */
+void ?{}(double& a, int b) { a = b; }
+double ?=?(double& n, zero_t) { return n = 0L; }
+// long double
+void ?{}(long double& a, int b) { a = b; }
+long double ?=?(long double& n, zero_t) { return n = 0L; }
+}
+
+trait dottable(otype V, otype T) {
+    T dot(V, V);
+};
+
+static inline {
+
+forall(otype T | sqrt(T), otype V | dottable(V, T))
+T length(V v) {
+   return sqrt(dot(v, v));
+}
+
+forall(otype T, otype V | dottable(V, T))
+T length_squared(V v) {
+   return dot(v, v);
+}
+
+forall(otype T, otype V | { T length(V); } | subtract(V))
+T distance(V v1, V v2) {
+    return length(v1 - v2);
+}
+
+
+forall(otype T, otype V | { T length(V); V ?/?(V, T); })
+V normalize(V v) {
+    return v / length(v);
+}
+
+// Project vector u onto vector v
+forall(otype T, otype V | dottable(V, T) | { V normalize(V); V ?*?(V, T); })
+V project(V u, V v) {
+    V v_norm = normalize(v);
+    return v_norm * dot(u, v_norm);
+}
+
+// Reflect incident vector v with respect to surface with normal n
+forall(otype T | fromint(T), otype V | { V project(V, V); V ?*?(T, V); V ?-?(V,V); })
+V reflect(V v, V n) {
+    return v - (T){2} * project(v, n);
+}
+
+// Refract incident vector v with respect to surface with normal n
+// eta is the ratio of indices of refraction between starting material and
+// entering material (i.e., from air to water, eta = 1/1.33)
+// v and n must already be normalized
+forall(otype T | fromint(T) | subtract(T) | multiply(T) | add(T) | lessthan(T) | sqrt(T),
+       otype V | dottable(V, T) | { V ?*?(T, V); V ?-?(V,V); void ?{}(V&, zero_t); })
+V refract(V v, V n, T eta) {
+    T dotValue = dot(n, v);
+    T k = (T){1} - eta * eta * ((T){1} - dotValue * dotValue);
+    if (k < (T){0}) {
+        return 0;
+    }
+    return eta * v - (eta * dotValue + sqrt(k)) * n;
+}
+
+// Given a perturbed normal and a geometric normal,
+// flip the perturbed normal if the geometric normal is pointing away
+// from the observer.
+// n is the perturbed vector that we want to align
+// i is the incident vector
+// ng is the geometric normal of the surface
+/* forall(| add(T) | multiply(T) | lessthan(T) | fromint(T) | subtract(T)) */
+forall(otype T | lessthan(T) | zeroinit(T), otype V | dottable(V, T) | negate(V))
+V faceforward(V n, V i, V ng) {
+    return dot(ng, i) < (T){0} ? n : -n;
+}
+
+} // inline
Index: libcfa/src/vec/vec2.hfa
===================================================================
--- libcfa/src/vec/vec2.hfa	(revision 24443247fe0f583ea2d0e1f705d9edf17e22d31f)
+++ libcfa/src/vec/vec2.hfa	(revision 3376ec9679544e4a3573e3f6ba122456afec811f)
@@ -1,53 +1,7 @@
 #pragma once
-#include <math.hfa>
-#include <iostream.hfa>
 
-trait zeroassn(otype T) {
-    T ?=?(T&, zero_t);
-};
-trait fromint(otype T) {
-    void ?{}(T&, int);
-};
-trait zero_assign(otype T) {
-    T ?=?(T&, zero_t);
-};
-trait subtract(otype T) {
-    T ?-?(T, T);
-    T -?(T);
-};
-trait add(otype T) {
-    T ?+?(T, T);
-};
-trait multiply(otype T) {
-    T ?*?(T, T);
-};
-trait divide(otype T) {
-    T ?/?(T, T);
-};
-trait lessthan(otype T) {
-    int ?<?(T, T);
-};
-trait equality(otype T) {
-    int ?==?(T, T);
-};
-trait sqrt(otype T) {
-    T sqrt(T);
-};
+#include "vec.hfa"
 
-static inline {
-// int
-int ?=?(int& n, zero_t) { return n = 0.f; }
-/* float */
-void ?{}(float& a, int b) { a = b; }
-float ?=?(float& n, zero_t) { return n = 0.f; }
-/* double */
-void ?{}(double& a, int b) { a = b; }
-double ?=?(double& n, zero_t) { return n = 0L; }
-// long double
-void ?{}(long double& a, int b) { a = b; }
-long double ?=?(long double& n, zero_t) { return n = 0L; }
-}
-
-forall(otype T) {
+forall (otype T) {
     struct vec2 {
         T x, y;
@@ -55,8 +9,7 @@
 }
 
-forall(otype T) {
+
+forall (otype T) {
     static inline {
-
-    // Constructors
 
     void ?{}(vec2(T)& v, T x, T y) {
@@ -98,4 +51,7 @@
         return u;
     }
+    }
+
+    forall(| negate(T)) {
     vec2(T) -?(vec2(T)& v) with (v) {
         return [-x, -y];
@@ -149,4 +105,5 @@
     }
 
+    // Geometric functions
     forall(| add(T) | multiply(T))
     T dot(vec2(T) u, vec2(T) v) {
@@ -154,71 +111,4 @@
     }
 
-    forall(| sqrt(T) | add(T) | multiply(T))
-    T length(vec2(T) v) {
-       return sqrt(dot(v, v));
-    }
-
-    forall(| add(T) | multiply(T))
-    T length_squared(vec2(T) v) {
-       return dot(v, v);
-    }
-
-    forall(| subtract(T) | sqrt(T) | add(T) | multiply(T))
-    T distance(vec2(T) v1, vec2(T) v2) {
-        return length(v1 - v2);
-    }
-
-    forall(| sqrt(T) | divide(T) | add(T) | multiply(T))
-    vec2(T) normalize(vec2(T) v) {
-        return v / sqrt(dot(v, v));
-    }
-
-    // Project vector u onto vector v
-    forall(| sqrt(T) | divide(T) | add(T) | multiply(T))
-    vec2(T) project(vec2(T) u, vec2(T) v) {
-        vec2(T) v_norm = normalize(v);
-        return v_norm * dot(u, v_norm);
-    }
-
-    // Reflect incident vector v with respect to surface with normal n
-    forall(| sqrt(T) | divide(T) | add(T) | multiply(T) | subtract(T) | fromint(T))
-    vec2(T) reflect(vec2(T) v, vec2(T) n) {
-        return v - (T){2} * project(v, n);
-    }
-
-    // Refract incident vector v with respect to surface with normal n
-    // eta is the ratio of indices of refraction between starting material and
-    // entering material (i.e., from air to water, eta = 1/1.33)
-    // v and n must already be normalized
-    forall(| sqrt(T) | add(T) | multiply(T) | subtract(T) | fromint(T) | lessthan(T) | zeroassn(T))
-    vec2(T) refract(vec2(T) v, vec2(T) n, T eta) {
-        T dotValue = dot(n, v);
-        T k = (T){1} - eta * eta * ((T){1} - dotValue * dotValue);
-        if (k < (T){0}) {
-            return 0;
-        }
-        return eta * v - (eta * dotValue + sqrt(k)) * n;
-    }
-
-    // Given a perturbed normal and a geometric normal,
-    // flip the perturbed normal if the geometric normal is pointing away
-    // from the observer.
-    // n is the perturbed vector that we want to align
-    // i is the incident vector
-    // ng is the geometric normal of the surface
-    forall(| add(T) | multiply(T) | lessthan(T) | fromint(T) | subtract(T))
-    vec2(T) faceforward(vec2(T) n, vec2(T) i, vec2(T) ng) {
-        return dot(ng, i) < (T){0} ? n : -n;
-    }
-
-    }
+    } // static inline
 }
-
-forall(dtype ostype, otype T | writeable(T, ostype)) {
-    ostype & ?|?( ostype & os, vec2(T) v) with (v) {
-        return os | '<' | x | ',' | y | '>';
-    }
-    void ?|?( ostype & os, vec2(T) v ) with (v) {
-        (ostype &)(os | v); ends(os);
-    }
-}
Index: bcfa/src/vec/vec2_f.hfa
===================================================================
--- libcfa/src/vec/vec2_f.hfa	(revision 24443247fe0f583ea2d0e1f705d9edf17e22d31f)
+++ 	(revision )
@@ -1,155 +1,0 @@
-#pragma once
-#include <math.hfa>
-/* #include <iostream.hfa> */
-
-struct vec2_f {
-    float x, y;
-};
-
-static inline {
-
-// Constructors
-
-void ?{}(vec2_f& v, float x, float y) {
-    v.[x, y] = [x, y];
-}
-void ?{}(vec2_f& vec, zero_t) with (vec) {
-    x = y = 0;
-}
-void ?{}(vec2_f& vec, float val) with (vec) {
-    x = y = val;
-}
-void ?{}(vec2_f& vec, const vec2_f& other) with (vec) {
-    [x,y] = other.[x,y];
-}
-
-// Assignment
-void ?=?(vec2_f& vec, const vec2_f& other) with (vec) {
-    [x,y] = other.[x,y];
-}
-void ?=?(vec2_f& vec, zero_t) with (vec) {
-    [x,y] = [0,0];
-}
-
-// Primitive mathematical operations
-
-// Subtraction
-vec2_f ?-?(const vec2_f& u, const vec2_f& v) {
-    return [u.x - v.x, u.y - v.y];
-}
-vec2_f& ?-=?(vec2_f& u, const vec2_f& v) {
-    u = u - v;
-    return u;
-}
-vec2_f -?(const vec2_f& v) with (v) {
-    return [-x, -y];
-}
-
-// Addition
-vec2_f ?+?(const vec2_f& u, const vec2_f& v) {
-    return [u.x + v.x, u.y + v.y];
-}
-vec2_f& ?+=?(vec2_f& u, const vec2_f& v) {
-    u = u + v;
-    return u;
-}
-
-// Scalar Multiplication
-vec2_f ?*?(const vec2_f& v, float scalar) with (v) {
-    return [x * scalar, y * scalar];
-}
-vec2_f ?*?(float scalar, const vec2_f& v) {
-    return v * scalar;
-}
-vec2_f& ?*=?(vec2_f& v, float scalar) with (v) {
-    v = v * scalar;
-    return v;
-}
-
-
-// Scalar Division
-vec2_f ?/?(const vec2_f& v, float scalar) with (v) {
-    return [x / scalar, y / scalar];
-}
-vec2_f& ?/=?(vec2_f& v, float scalar) with (v) {
-    v = v / scalar;
-    return v;
-}
-
-// Relational Operators
-bool ?==?(const vec2_f& u, const vec2_f& v) with (u) {
-    return x == v.x && y == v.y;
-}
-bool ?!=?(const vec2_f& u, const vec2_f& v) {
-    return !(u == v);
-}
-
-/* // Printing the vector (ostream) */
-/* forall( dtype ostype | ostream( ostype ) ) { */
-/*     ostype & ?|?( ostype & os, const vec2_f& v) with (v) { */
-/*         if ( sepPrt( os ) ) fmt( os, "%s", sepGetCur( os ) ); */
-/*         fmt( os, "<%g,%g>", x, y); */
-/*         return os; */
-/*     } */
-/*     void ?|?( ostype & os, const vec2_f& v ) { */
-/*         (ostype &)(os | v); ends( os ); */
-/*     } */
-/* } */
-
-//---------------------- Geometric Functions ----------------------
-// These functions implement the Geometric Functions section of GLSL for 2D vectors
-
-float dot(const vec2_f& u, const vec2_f& v) {
-    return u.x * v.x + u.y * v.y;
-}
-
-float length(const vec2_f& v) {
-   return sqrt(dot(v, v));
-}
-
-float length_squared(const vec2_f& v) {
-   return dot(v, v);
-}
-
-float distance(const vec2_f& v1, const vec2_f& v2) {
-    return length(v1 - v2);
-}
-
-vec2_f normalize(const vec2_f& v) {
-    return v / sqrt(dot(v, v));
-}
-
-// Project vector u onto vector v
-vec2_f project(const vec2_f& u, const vec2_f& v) {
-    vec2_f v_norm = normalize(v);
-    return v_norm * dot(u, v_norm);
-}
-
-// Reflect incident vector v with respect to surface with normal n
-vec2_f reflect(const vec2_f& v, const vec2_f& n) {
-    return v - 2 * project(v, n);
-}
-
-// Refract incident vector v with respect to surface with normal n
-// eta is the ratio of indices of refraction between starting material and
-// entering material (i.e., from air to water, eta = 1/1.33)
-vec2_f refract(const vec2_f& v, const vec2_f& n, float eta) {
-    float dotValue = dot(n, v);
-    float k = 1 - eta \ 2 * (1 - dotValue \ 2);
-    if (k < 0) {
-        return 0;
-    }
-    return eta * v - (eta * dotValue + sqrt(k)) * n;
-}
-
-// Given a perturbed normal and a geometric normal,
-// flip the perturbed normal if the geometric normal is pointing away
-// from the observer.
-// n is the perturbed vector that we want to align
-// i is the incident vector
-// ng is the geometric normal of the surface
-vec2_f faceforward(const vec2_f& n, const vec2_f& i, const vec2_f& ng) {
-    return dot(ng, i) < 0 ? n : -n;
-}
-
-}
Index: bcfa/src/vec/vec3.hfa
===================================================================
--- libcfa/src/vec/vec3.hfa	(revision 24443247fe0f583ea2d0e1f705d9edf17e22d31f)
+++ 	(revision )
@@ -1,187 +1,0 @@
-#pragma once
-#include <math.hfa>
-#include <iostream.hfa>
-
-trait vec3_t(otype T) {
-    void ?{}(T&, int);
-    T ?=?(T&, zero_t);
-    T ?-?(T, T);
-    T -?(T);
-    T ?+?(T, T);
-    T ?*?(T, T);
-    T ?/?(T, T);
-    int ?==?(T, T);
-    int ?<?(T, T);
-    T sqrt(T);
-};
-
-static inline {
-// int
-int ?=?(int& n, zero_t) { return n = 0.f; }
-int sqrt(int a) { return sqrt((float)a); }
-/* float */
-void ?{}(float& a, int b) { a = b; }
-float ?=?(float& n, zero_t) { return n = 0.f; }
-/* double */
-void ?{}(double& a, int b) { a = b; }
-double ?=?(double& n, zero_t) { return n = 0L; }
-// long double
-void ?{}(long double& a, int b) { a = b; }
-long double ?=?(long double& n, zero_t) { return n = 0L; }
-}
-
-forall(otype T | vec3_t(T)) {
-    struct vec3 {
-        T x, y, z;
-    };
-}
-
-/* static inline { */
-forall(otype T | vec3_t(T)) {
-    static inline {
-
-    // Constructors
-
-    void ?{}(vec3(T)& v, T x, T y, T z) {
-        v.[x, y, z] = [x, y, z];
-    }
-    void ?{}(vec3(T)& vec, zero_t) with (vec) {
-        x = y = z = 0;
-    }
-    void ?{}(vec3(T)& vec, T val) with (vec) {
-        x = y = z = val;
-    }
-    void ?{}(vec3(T)& vec, vec3(T) other) with (vec) {
-        [x,y,z] = other.[x,y,z];
-    }
-
-    // Assignment
-    void ?=?(vec3(T)& vec, vec3(T) other) with (vec) {
-        [x,y,z] = other.[x,y,z];
-    }
-    void ?=?(vec3(T)& vec, zero_t) with (vec) {
-        x = y = z = 0;
-    }
-
-    // Primitive mathematical operations
-
-    // Subtraction
-    vec3(T) ?-?(vec3(T) u, vec3(T) v) { // TODO( can't make this const ref )
-        return [u.x - v.x, u.y - v.y, u.z - v.z];
-    }
-    vec3(T)& ?-=?(vec3(T)& u, vec3(T) v) {
-        u = u - v;
-        return u;
-    }
-    vec3(T) -?(vec3(T)& v) with (v) {
-        return [-x, -y, -z];
-    }
-
-    // Addition
-    vec3(T) ?+?(vec3(T) u, vec3(T) v) { // TODO( can't make this const ref )
-        return [u.x + v.x, u.y + v.y, u.z + v.z];
-    }
-    vec3(T)& ?+=?(vec3(T)& u, vec3(T) v) {
-        u = u + v;
-        return u;
-    }
-
-    // Scalar Multiplication
-    vec3(T) ?*?(vec3(T) v, T scalar) with (v) { // TODO (can't make this const ref)
-        return [x * scalar, y * scalar, z * scalar];
-    }
-    vec3(T) ?*?(T scalar, vec3(T) v) { // TODO (can't make this const ref)
-        return v * scalar;
-    }
-    vec3(T)& ?*=?(vec3(T)& v, T scalar) {
-        v = v * scalar;
-        return v;
-    }
-
-
-    // Scalar Division
-    vec3(T) ?/?(vec3(T) v, T scalar) with (v) {
-        return [x / scalar, y / scalar, z / scalar];
-    }
-    vec3(T)& ?/=?(vec3(T)& v, T scalar) with (v) {
-        v = v / scalar;
-        return v;
-    }
-    // Relational Operators
-    bool ?==?(vec3(T) u, vec3(T) v) with (u) {
-        return x == v.x && y == v.y && z == v.z;
-    }
-    bool ?!=?(vec3(T) u, vec3(T) v) {
-        return !(u == v);
-    }
-
-    T dot(vec3(T) u, vec3(T) v) {
-        return u.x * v.x + u.y * v.y + u.z * v.z;
-    }
-
-    T length(vec3(T) v) {
-       return sqrt(dot(v, v));
-    }
-
-    T length_squared(vec3(T) v) {
-       return dot(v, v);
-    }
-
-    T distance(vec3(T) v1, vec3(T) v2) {
-        return length(v1 - v2);
-    }
-
-    vec3(T) normalize(vec3(T) v) {
-        return v / sqrt(dot(v, v));
-    }
-
-    // Project vector u onto vector v
-    vec3(T) project(vec3(T) u, vec3(T) v) {
-        vec3(T) v_norm = normalize(v);
-        return v_norm * dot(u, v_norm);
-    }
-
-    // Reflect incident vector v with respect to surface with normal n
-    vec3(T) reflect(vec3(T) v, vec3(T) n) {
-        return v - (T){2} * project(v, n);
-    }
-
-    // Refract incident vector v with respect to surface with normal n
-    // eta is the ratio of indices of refraction between starting material and
-    // entering material (i.e., from air to water, eta = 1/1.33)
-    // v and n must already be normalized
-    vec3(T) refract(vec3(T) v, vec3(T) n, T eta) {
-        T dotValue = dot(n, v);
-        T k = (T){1} - eta * eta * ((T){1} - dotValue * dotValue);
-        if (k < (T){0}) {
-            return 0;
-        }
-        return eta * v - (eta * dotValue + sqrt(k)) * n;
-    }
-
-    // Given a perturbed normal and a geometric normal,
-    // flip the perturbed normal if the geometric normal is pointing away
-    // from the observer.
-    // n is the perturbed vector that we want to align
-    // i is the incident vector
-    // ng is the geometric normal of the surface
-    vec3(T) faceforward(vec3(T) n, vec3(T) i, vec3(T) ng) {
-        return dot(ng, i) < (T){0} ? n : -n;
-    }
-
-    vec3(T) cross(vec3(T) u, vec3(T) v) {
-        return (vec3(T)){ u.y * v.z - v.y * u.z,
-                          u.z * v.x - v.z * u.x,
-                          u.x * v.y - v.x * u.y };
-    }
-    }
-}
-
-forall(dtype ostype, otype T | writeable(T, ostype) | vec3_t(T)) {
-    ostype & ?|?( ostype & os, vec3(T) v) with (v) {
-        return os | '<' | x | ',' | y | ',' | z | '>';
-    }
-    void ?|?( ostype & os, vec3(T) v ) with (v) {
-        (ostype &)(os | v); ends(os);
-    }
-}