summaryrefslogtreecommitdiffstats
path: root/libdimension/dimension/math
diff options
context:
space:
mode:
Diffstat (limited to 'libdimension/dimension/math')
-rw-r--r--libdimension/dimension/math/aabb.h25
-rw-r--r--libdimension/dimension/math/matrix.h210
-rw-r--r--libdimension/dimension/math/vector.h104
3 files changed, 213 insertions, 126 deletions
diff --git a/libdimension/dimension/math/aabb.h b/libdimension/dimension/math/aabb.h
index 14cc575..4eb7870 100644
--- a/libdimension/dimension/math/aabb.h
+++ b/libdimension/dimension/math/aabb.h
@@ -57,8 +57,8 @@ DMNSN_INLINE dmnsn_aabb
dmnsn_zero_aabb(void)
{
dmnsn_aabb box = {
- { DMNSN_INFINITY, DMNSN_INFINITY, DMNSN_INFINITY },
- { -DMNSN_INFINITY, -DMNSN_INFINITY, -DMNSN_INFINITY }
+ { { DMNSN_INFINITY, DMNSN_INFINITY, DMNSN_INFINITY } },
+ { { -DMNSN_INFINITY, -DMNSN_INFINITY, -DMNSN_INFINITY } }
};
return box;
}
@@ -68,8 +68,8 @@ DMNSN_INLINE dmnsn_aabb
dmnsn_infinite_aabb(void)
{
dmnsn_aabb box = {
- { -DMNSN_INFINITY, -DMNSN_INFINITY, -DMNSN_INFINITY },
- { DMNSN_INFINITY, DMNSN_INFINITY, DMNSN_INFINITY }
+ { { -DMNSN_INFINITY, -DMNSN_INFINITY, -DMNSN_INFINITY } },
+ { { DMNSN_INFINITY, DMNSN_INFINITY, DMNSN_INFINITY } }
};
return box;
}
@@ -94,15 +94,26 @@ dmnsn_symmetric_aabb(dmnsn_vector r)
DMNSN_INLINE bool
dmnsn_aabb_contains(dmnsn_aabb box, dmnsn_vector p)
{
- return (p.x >= box.min.x && p.y >= box.min.y && p.z >= box.min.z)
- && (p.x <= box.max.x && p.y <= box.max.y && p.z <= box.max.z);
+ for (unsigned int i = 0; i < 3; ++i) {
+ if (p.n[i] < box.min.n[i]) {
+ return false;
+ }
+ }
+
+ for (unsigned int i = 0; i < 3; ++i) {
+ if (p.n[i] > box.max.n[i]) {
+ return false;
+ }
+ }
+
+ return true;
}
/** Return whether a bounding box is infinite. */
DMNSN_INLINE bool
dmnsn_aabb_is_infinite(dmnsn_aabb box)
{
- return box.min.x == -DMNSN_INFINITY;
+ return box.min.n[0] == -DMNSN_INFINITY;
}
/**
diff --git a/libdimension/dimension/math/matrix.h b/libdimension/dimension/math/matrix.h
index 7471bf5..e67121e 100644
--- a/libdimension/dimension/math/matrix.h
+++ b/libdimension/dimension/math/matrix.h
@@ -39,152 +39,218 @@ typedef struct dmnsn_matrix {
"[%g\t%g\t%g\t%g]\n" \
"[%g\t%g\t%g\t%g]"
/** The appropriate arguements to printf() a matrix. */
-#define DMNSN_MATRIX_PRINTF(m) \
- (m).n[0][0], (m).n[0][1], (m).n[0][2], (m).n[0][3], \
- (m).n[1][0], (m).n[1][1], (m).n[1][2], (m).n[1][3], \
- (m).n[2][0], (m).n[2][1], (m).n[2][2], (m).n[2][3], \
+#define DMNSN_MATRIX_PRINTF(M) \
+ (M).n[0][0], (M).n[0][1], (M).n[0][2], (M).n[0][3], \
+ (M).n[1][0], (M).n[1][1], (M).n[1][2], (M).n[1][3], \
+ (M).n[2][0], (M).n[2][1], (M).n[2][2], (M).n[2][3], \
0.0, 0.0, 0.0, 1.0
-/** Construct a new transformation matrix. */
+/** Create a transformation matrix. */
DMNSN_INLINE dmnsn_matrix
-dmnsn_new_matrix(double a0, double a1, double a2, double a3,
- double b0, double b1, double b2, double b3,
- double c0, double c1, double c2, double c3)
+dmnsn_new_matrix(double a0, double b0, double c0, double d0,
+ double a1, double b1, double c1, double d1,
+ double a2, double b2, double c2, double d2)
{
- dmnsn_matrix m = { { { a0, a1, a2, a3 },
- { b0, b1, b2, b3 },
- { c0, c1, c2, c3 } } };
- return m;
+ dmnsn_matrix M = {
+ {
+ { a0, b0, c0, d0 },
+ { a1, b1, c1, d1 },
+ { a2, b2, c2, d2 }
+ }
+ };
+ return M;
}
-/** Construct a new transformation matrix from column vectors. */
+/** Create a transformation matrix from column vectors. */
DMNSN_INLINE dmnsn_matrix
-dmnsn_new_matrix4(dmnsn_vector a, dmnsn_vector b, dmnsn_vector c,
- dmnsn_vector d)
+dmnsn_new_matrix4(dmnsn_vector a, dmnsn_vector b, dmnsn_vector c, dmnsn_vector d)
{
- dmnsn_matrix m = { { { a.x, b.x, c.x, d.x },
- { a.y, b.y, c.y, d.y },
- { a.z, b.z, c.z, d.z } } };
- return m;
+ dmnsn_matrix M;
+
+ unsigned int i;
+ for (i = 0; i < 3; ++i) {
+ M.n[i][0] = a.n[i];
+ }
+ for (i = 0; i < 3; ++i) {
+ M.n[i][1] = b.n[i];
+ }
+ for (i = 0; i < 3; ++i) {
+ M.n[i][2] = c.n[i];
+ }
+ for (i = 0; i < 3; ++i) {
+ M.n[i][3] = d.n[i];
+ }
+
+ return M;
}
/** Extract column vectors from a matrix. */
DMNSN_INLINE dmnsn_vector
dmnsn_matrix_column(dmnsn_matrix M, unsigned int i)
{
- return dmnsn_new_vector(M.n[0][i], M.n[1][i], M.n[2][i]);
+ dmnsn_vector v;
+ unsigned int j;
+ for (j = 0; j < 3; ++j) {
+ v.n[j] = M.n[j][i];
+ }
+ return v;
}
/** Return the identity matrix. */
-dmnsn_matrix dmnsn_identity_matrix(void);
+DMNSN_INLINE dmnsn_matrix
+dmnsn_identity_matrix(void)
+{
+ return dmnsn_new_matrix(
+ 1.0, 0.0, 0.0, 0.0,
+ 0.0, 1.0, 0.0, 0.0,
+ 0.0, 0.0, 1.0, 0.0
+ );
+}
/**
- * A scale transformation.
- * @param[in] s A vector with components representing the scaling factor in
- * each axis.
- * @return The transformation matrix.
+ * Return a scale transformation.
+ * @param[in] s A vector with components representing the scaling factor in
+ * each axis.
*/
-dmnsn_matrix dmnsn_scale_matrix(dmnsn_vector s);
+DMNSN_INLINE dmnsn_matrix
+dmnsn_scale_matrix(dmnsn_vector s)
+{
+ return dmnsn_new_matrix(
+ s.n[0], 0.0, 0.0, 0.0,
+ 0.0, s.n[1], 0.0, 0.0,
+ 0.0, 0.0, s.n[2], 0.0
+ );
+}
+
/**
- * A translation.
- * @param[in] d The vector to translate by.
- * @return The transformation matrix.
+ * Set \p M to a translation matrix.
+ * @param[in] d The vector to translate by.
*/
-dmnsn_matrix dmnsn_translation_matrix(dmnsn_vector d);
+DMNSN_INLINE dmnsn_matrix
+dmnsn_translation_matrix(dmnsn_vector d)
+{
+ return dmnsn_new_matrix(
+ 1.0, 0.0, 0.0, d.n[0],
+ 0.0, 1.0, 0.0, d.n[1],
+ 0.0, 0.0, 1.0, d.n[2]
+ );
+}
+
/**
- * A left-handed rotation.
- * @param[in] theta A vector representing an axis and angle.
- * @f$ axis = \vec{\theta}/|\vec{\theta}| @f$,
- * @f$ angle = |\vec{\theta}| @f$
- * @return The transformation matrix.
+ * Return a rotation matrix.
+ * @param[in] theta A vector representing an axis and angle.
+ * @f$ axis = \vec{\theta}/|\vec{\theta}| @f$,
+ * @f$ angle = |\vec{\theta}| @f$
*/
dmnsn_matrix dmnsn_rotation_matrix(dmnsn_vector theta);
/**
- * An alignment matrix.
+ * Return an alignment matrix.
* @param[in] from The initial vector.
* @param[in] to The desired direction.
* @param[in] axis1 The first axis about which to rotate.
* @param[in] axis2 The second axis about which to rotate.
- * @return A transformation matrix that will rotate \p from to \p to.
*/
-dmnsn_matrix dmnsn_alignment_matrix(dmnsn_vector from, dmnsn_vector to,
- dmnsn_vector axis1, dmnsn_vector axis2);
+dmnsn_matrix dmnsn_alignment_matrix(dmnsn_vector from, dmnsn_vector to, dmnsn_vector axis1, dmnsn_vector axis2);
/** Invert a matrix. */
-dmnsn_matrix dmnsn_matrix_inverse(dmnsn_matrix A);
+dmnsn_matrix dmnsn_matrix_inverse(dmnsn_matrix M);
/** Multiply two matricies. */
-dmnsn_matrix dmnsn_matrix_mul(dmnsn_matrix lhs, dmnsn_matrix rhs);
+DMNSN_INLINE dmnsn_matrix dmnsn_matrix_mul(dmnsn_matrix lhs, dmnsn_matrix rhs)
+{
+ dmnsn_matrix M;
+
+ unsigned int i, j, k;
+ for (i = 0; i < 3; ++i) {
+ for (j = 0; j < 3; ++j) {
+ M.n[i][j] = 0.0;
+ for (k = 0; k < 3; ++k) {
+ M.n[i][j] += lhs.n[i][k]*rhs.n[k][j];
+ }
+ }
+
+ M.n[i][3] = lhs.n[i][3];
+ for (k = 0; k < 3; ++k) {
+ M.n[i][3] += lhs.n[i][k]*rhs.n[k][3];
+ }
+ }
+
+ return M;
+}
/** Transform a point by a matrix. */
DMNSN_INLINE dmnsn_vector
-dmnsn_transform_point(dmnsn_matrix T, dmnsn_vector v)
+dmnsn_transform_point(dmnsn_matrix M, dmnsn_vector v)
{
- /* 9 multiplications, 9 additions */
dmnsn_vector r;
- r.x = T.n[0][0]*v.x + T.n[0][1]*v.y + T.n[0][2]*v.z + T.n[0][3];
- r.y = T.n[1][0]*v.x + T.n[1][1]*v.y + T.n[1][2]*v.z + T.n[1][3];
- r.z = T.n[2][0]*v.x + T.n[2][1]*v.y + T.n[2][2]*v.z + T.n[2][3];
+ unsigned int i, j;
+ for (i = 0; i < 3; ++i) {
+ r.n[i] = M.n[i][3];
+ for (j = 0; j < 3; ++j) {
+ r.n[i] += M.n[i][j]*v.n[j];
+ }
+ }
return r;
}
/** Transform a direction by a matrix. */
DMNSN_INLINE dmnsn_vector
-dmnsn_transform_direction(dmnsn_matrix T, dmnsn_vector v)
+dmnsn_transform_direction(dmnsn_matrix M, dmnsn_vector v)
{
- /* 9 multiplications, 6 additions */
dmnsn_vector r;
- r.x = T.n[0][0]*v.x + T.n[0][1]*v.y + T.n[0][2]*v.z;
- r.y = T.n[1][0]*v.x + T.n[1][1]*v.y + T.n[1][2]*v.z;
- r.z = T.n[2][0]*v.x + T.n[2][1]*v.y + T.n[2][2]*v.z;
+ unsigned int i, j;
+ for (i = 0; i < 3; ++i) {
+ r.n[i] = 0.0;
+ for (j = 0; j < 3; ++j) {
+ r.n[i] += M.n[i][j]*v.n[j];
+ }
+ }
return r;
}
/**
* Transform a pseudovector by a matrix.
- * @param[in] Tinv The inverse of the transformation matrix.
- * @param[in] v The pseudovector to transform
+ * @param[in] Minv The inverse of the transformation matrix.
+ * @param[in] v The pseudovector to transform.
* @return The transformed pseudovector.
*/
DMNSN_INLINE dmnsn_vector
-dmnsn_transform_normal(dmnsn_matrix Tinv, dmnsn_vector v)
+dmnsn_transform_normal(dmnsn_matrix Minv, dmnsn_vector v)
{
- /* Multiply by the transpose of the inverse
- (9 multiplications, 6 additions) */
+ /* Multiply by the transpose of the inverse */
dmnsn_vector r;
- r.x = Tinv.n[0][0]*v.x + Tinv.n[1][0]*v.y + Tinv.n[2][0]*v.z;
- r.y = Tinv.n[0][1]*v.x + Tinv.n[1][1]*v.y + Tinv.n[2][1]*v.z;
- r.z = Tinv.n[0][2]*v.x + Tinv.n[1][2]*v.y + Tinv.n[2][2]*v.z;
+ unsigned int i, j;
+ for (i = 0; i < 3; ++i) {
+ r.n[i] = 0.0;
+ for (j = 0; j < 3; ++j) {
+ r.n[i] += Minv.n[j][i]*v.n[j];
+ }
+ }
return r;
}
-/**
- * Transform a ray by a matrix.
- * \f$ n' = T(l.\vec{x_0} + l.\vec{n}) - T(l.\vec{x_0}) \f$,
- * \f$ \vec{x_0}' = T(l.\vec{x_0}) \f$
- */
+/** Transform a ray by a matrix. */
DMNSN_INLINE dmnsn_ray
-dmnsn_transform_ray(dmnsn_matrix T, dmnsn_ray l)
+dmnsn_transform_ray(dmnsn_matrix M, dmnsn_ray l)
{
- /* 18 multiplications, 15 additions */
dmnsn_ray ret;
- ret.x0 = dmnsn_transform_point(T, l.x0);
- ret.n = dmnsn_transform_direction(T, l.n);
+ ret.x0 = dmnsn_transform_point(M, l.x0);
+ ret.n = dmnsn_transform_direction(M, l.n);
return ret;
}
/** Transform a bounding box by a matrix. */
-dmnsn_aabb dmnsn_transform_aabb(dmnsn_matrix T, dmnsn_aabb box);
+dmnsn_aabb dmnsn_transform_aabb(dmnsn_matrix M, dmnsn_aabb box);
/** Return whether a matrix contains any NaN components. */
DMNSN_INLINE bool
-dmnsn_matrix_isnan(dmnsn_matrix m)
+dmnsn_matrix_isnan(dmnsn_matrix M)
{
- size_t i, j;
+ unsigned int i, j;
for (i = 0; i < 3; ++i) {
for (j = 0; j < 4; ++j) {
- if (dmnsn_isnan(m.n[i][j])) {
+ if (dmnsn_isnan(M.n[i][j])) {
return true;
}
}
diff --git a/libdimension/dimension/math/vector.h b/libdimension/dimension/math/vector.h
index 8eacee9..90b9a3d 100644
--- a/libdimension/dimension/math/vector.h
+++ b/libdimension/dimension/math/vector.h
@@ -32,26 +32,24 @@
/** A vector in 3 dimensions. */
typedef struct dmnsn_vector {
- double x; /**< The x component. */
- double y; /**< The y component. */
- double z; /**< The z component. */
+ double n[3]; /**< The components. */
} dmnsn_vector;
/** A standard format string for vectors. */
#define DMNSN_VECTOR_FORMAT "<%g, %g, %g>"
/** The appropriate arguements to printf() a vector. */
-#define DMNSN_VECTOR_PRINTF(v) (v).x, (v).y, (v).z
+#define DMNSN_VECTOR_PRINTF(v) (v).n[0], (v).n[1], (v).n[2]
/* Constants */
/** The zero vector. */
-static const dmnsn_vector dmnsn_zero = { 0.0, 0.0, 0.0 };
+static const dmnsn_vector dmnsn_zero = { { 0.0, 0.0, 0.0 } };
/** The x vector. */
-static const dmnsn_vector dmnsn_x = { 1.0, 0.0, 0.0 };
+static const dmnsn_vector dmnsn_x = { { 1.0, 0.0, 0.0 } };
/** The y vector. */
-static const dmnsn_vector dmnsn_y = { 0.0, 1.0, 0.0 };
+static const dmnsn_vector dmnsn_y = { { 0.0, 1.0, 0.0 } };
/** The z vector. */
-static const dmnsn_vector dmnsn_z = { 0.0, 0.0, 1.0 };
+static const dmnsn_vector dmnsn_z = { { 0.0, 0.0, 1.0 } };
/* Shorthand for vector construction */
@@ -59,7 +57,7 @@ static const dmnsn_vector dmnsn_z = { 0.0, 0.0, 1.0 };
DMNSN_INLINE dmnsn_vector
dmnsn_new_vector(double x, double y, double z)
{
- dmnsn_vector v = { x, y, z };
+ dmnsn_vector v = { { x, y, z } };
return v;
}
@@ -69,8 +67,11 @@ dmnsn_new_vector(double x, double y, double z)
DMNSN_INLINE dmnsn_vector
dmnsn_vector_negate(dmnsn_vector rhs)
{
- /* 3 negations */
- dmnsn_vector v = { -rhs.x, -rhs.y, -rhs.z };
+ dmnsn_vector v;
+ unsigned int i;
+ for (i = 0; i < 3; ++i) {
+ v.n[i] = -rhs.n[i];
+ }
return v;
}
@@ -78,8 +79,11 @@ dmnsn_vector_negate(dmnsn_vector rhs)
DMNSN_INLINE dmnsn_vector
dmnsn_vector_add(dmnsn_vector lhs, dmnsn_vector rhs)
{
- /* 3 additions */
- dmnsn_vector v = { lhs.x + rhs.x, lhs.y + rhs.y, lhs.z + rhs.z };
+ dmnsn_vector v;
+ unsigned int i;
+ for (i = 0; i < 3; ++i) {
+ v.n[i] = lhs.n[i] + rhs.n[i];
+ }
return v;
}
@@ -87,8 +91,11 @@ dmnsn_vector_add(dmnsn_vector lhs, dmnsn_vector rhs)
DMNSN_INLINE dmnsn_vector
dmnsn_vector_sub(dmnsn_vector lhs, dmnsn_vector rhs)
{
- /* 3 additions */
- dmnsn_vector v = { lhs.x - rhs.x, lhs.y - rhs.y, lhs.z - rhs.z };
+ dmnsn_vector v;
+ unsigned int i;
+ for (i = 0; i < 3; ++i) {
+ v.n[i] = lhs.n[i] - rhs.n[i];
+ }
return v;
}
@@ -96,17 +103,11 @@ dmnsn_vector_sub(dmnsn_vector lhs, dmnsn_vector rhs)
DMNSN_INLINE dmnsn_vector
dmnsn_vector_mul(double lhs, dmnsn_vector rhs)
{
- /* 3 multiplications */
- dmnsn_vector v = { lhs*rhs.x, lhs*rhs.y, lhs*rhs.z };
- return v;
-}
-
-/** Divide a vector by a scalar. */
-DMNSN_INLINE dmnsn_vector
-dmnsn_vector_div(dmnsn_vector lhs, double rhs)
-{
- /* 3 divisions */
- dmnsn_vector v = { lhs.x/rhs, lhs.y/rhs, lhs.z/rhs };
+ dmnsn_vector v;
+ unsigned int i;
+ for (i = 0; i < 3; ++i) {
+ v.n[i] = lhs*rhs.n[i];
+ }
return v;
}
@@ -114,18 +115,22 @@ dmnsn_vector_div(dmnsn_vector lhs, double rhs)
DMNSN_INLINE double
dmnsn_vector_dot(dmnsn_vector lhs, dmnsn_vector rhs)
{
- /* 3 multiplications, 2 additions */
- return lhs.x*rhs.x + lhs.y*rhs.y + lhs.z*rhs.z;
+ double result = 0.0;
+ unsigned int i;
+ for (i = 0; i < 3; ++i) {
+ result += lhs.n[i]*rhs.n[i];
+ }
+ return result;
}
/** Return the cross product of two vectors. */
DMNSN_INLINE dmnsn_vector
dmnsn_vector_cross(dmnsn_vector lhs, dmnsn_vector rhs)
{
- /* 6 multiplications, 3 additions */
- dmnsn_vector v = { lhs.y*rhs.z - lhs.z*rhs.y,
- lhs.z*rhs.x - lhs.x*rhs.z,
- lhs.x*rhs.y - lhs.y*rhs.x };
+ dmnsn_vector v;
+ v.n[0] = lhs.n[1]*rhs.n[2] - lhs.n[2]*rhs.n[1];
+ v.n[1] = lhs.n[2]*rhs.n[0] - lhs.n[0]*rhs.n[2];
+ v.n[2] = lhs.n[0]*rhs.n[1] - lhs.n[1]*rhs.n[0];
return v;
}
@@ -133,7 +138,6 @@ dmnsn_vector_cross(dmnsn_vector lhs, dmnsn_vector rhs)
DMNSN_INLINE dmnsn_vector
dmnsn_vector_proj(dmnsn_vector u, dmnsn_vector d)
{
- /* 1 division, 9 multiplications, 4 additions */
return dmnsn_vector_mul(dmnsn_vector_dot(u, d)/dmnsn_vector_dot(d, d), d);
}
@@ -141,7 +145,6 @@ dmnsn_vector_proj(dmnsn_vector u, dmnsn_vector d)
DMNSN_INLINE double
dmnsn_vector_norm(dmnsn_vector n)
{
- /* 1 sqrt, 3 multiplications, 2 additions */
return sqrt(dmnsn_vector_dot(n, n));
}
@@ -149,35 +152,42 @@ dmnsn_vector_norm(dmnsn_vector n)
DMNSN_INLINE dmnsn_vector
dmnsn_vector_normalized(dmnsn_vector n)
{
- /* 1 sqrt, 3 divisions, 3 multiplications, 2 additions */
- return dmnsn_vector_div(n, dmnsn_vector_norm(n));
+ return dmnsn_vector_mul(1.0/dmnsn_vector_norm(n), n);
}
/** Return the component-wise minimum of two vectors. */
DMNSN_INLINE dmnsn_vector
dmnsn_vector_min(dmnsn_vector a, dmnsn_vector b)
{
- return dmnsn_new_vector(
- dmnsn_min(a.x, b.x),
- dmnsn_min(a.y, b.y),
- dmnsn_min(a.z, b.z)
- );
+ dmnsn_vector v;
+ unsigned int i;
+ for (i = 0; i < 3; ++i) {
+ v.n[i] = dmnsn_min(a.n[i], b.n[i]);
+ }
+ return v;
}
/** Return the component-wise maximum of two vectors. */
DMNSN_INLINE dmnsn_vector
dmnsn_vector_max(dmnsn_vector a, dmnsn_vector b)
{
- return dmnsn_new_vector(
- dmnsn_max(a.x, b.x),
- dmnsn_max(a.y, b.y),
- dmnsn_max(a.z, b.z)
- );
+ dmnsn_vector v;
+ unsigned int i;
+ for (i = 0; i < 3; ++i) {
+ v.n[i] = dmnsn_max(a.n[i], b.n[i]);
+ }
+ return v;
}
/** Return whether a vector contains any NaN components. */
DMNSN_INLINE bool
dmnsn_vector_isnan(dmnsn_vector v)
{
- return dmnsn_isnan(v.x) || dmnsn_isnan(v.y) || dmnsn_isnan(v.z);
+ unsigned int i;
+ for (i = 0; i < 3; ++i) {
+ if (dmnsn_isnan(v.n[i])) {
+ return true;
+ }
+ }
+ return false;
}