Import manifold code from fieldviz

4 files changed, 263 insertions, 0 deletions

diff --git a/include/fml/curve.h b/include/fml/curve.h
new file mode 100644
index 0000000..8a705cd
--- /dev/null
+++ b/include/fml/curve.h
@@ -0,0 +1,113 @@
+#ifndef CURVE_H
+#define CURVE_H
+
+#include <cmath>
+#include <iostream>
+
+#include "fml.h"
+
+#include "manifold.h"
+
+namespace fml {
+
+/* All curves inherit this class (which is empty because Manifold has
+ * everything we need). */
+    class Curve : public Manifold {
+    public:
+        const int dimension() const { return 1; }
+    };
+
+    class LineSegment : public Curve {
+    private:
+        vec3 a, b;
+    public:
+        LineSegment(vec3 a_, vec3 b_) : a(a_), b(b_) {};
+
+        vec3 integrate(vec3 (*integrand)(vec3 s, vec3 ds), scalar delta) const;
+
+        const char *name() const { return "LineSegment"; }
+
+        friend std::istream operator>>(std::istream &is, LineSegment &ls);
+    };
+
+    std::istream operator>>(std::istream &is, LineSegment &ls);
+
+    class Arc : public Curve {
+    private:
+        vec3 center;
+
+        /* these are relative to the center (direction will be determined
+         * by RHR of normal), and should be orthonormal */
+        vec3 radius, normal;
+
+        /* how many radians the arc extends for (can be greater than 2pi) */
+        scalar angle;
+    public:
+        Arc(vec3 c_, vec3 r_, vec3 n_, scalar th) : center(c_), radius(r_), normal(n_), angle(th) {};
+
+        vec3 integrate(vec3 (*integrand)(vec3 s, vec3 ds), scalar delta) const;
+
+        const char *name() const { return "Arc"; }
+
+        friend std::istream operator>>(std::istream &is, LineSegment &ls);
+    };
+
+    std::istream operator>>(std::istream &is, LineSegment &ls);
+
+    class Spiral : public Curve {
+    private:
+        vec3 origin;
+
+        /* these are relative to the center (direction will be determined
+         * by RHR of normal), and should be orthonormal */
+        vec3 radius, normal;
+
+        /* how many radians the arc extends for (can be greater than 2pi) */
+        scalar angle;
+
+        /* linear distance between turns (2pi) */
+        scalar pitch;
+    public:
+        Spiral(vec3 c_, vec3 r_, vec3 n_, scalar th, scalar p) : origin(c_), radius(r_), normal(n_), angle(th), pitch(p) {};
+
+        vec3 integrate(vec3 (*integrand)(vec3 s, vec3 ds), scalar delta) const;
+
+        const char *name() const { return "Solenoid"; }
+
+        friend std::istream operator>>(std::istream &is, LineSegment &ls);
+    };
+
+    std::istream operator>>(std::istream &is, LineSegment &ls);
+
+    class Toroid : public Curve {
+    private:
+        vec3 origin;
+
+        /* these are relative to the center (direction will be determined
+         * by RHR of normal), and should be orthonormal */
+        vec3 major_radius, major_normal;
+
+        /* "thickness" of toroid */
+        scalar minor_r;
+
+        /* how many radians (about the center) the toroid extends for
+         * (can be greater than 2pi) */
+        scalar major_angle;
+
+        /* central angle between successive turns (2pi rotation of small
+         * radius vector) */
+        scalar pitch;
+    public:
+        Toroid() {};
+        Toroid(vec3 o, vec3 maj_r, vec3 maj_n, scalar ang, scalar min_r, scalar p) : origin(o), major_radius(maj_r), major_normal(maj_n), major_angle(ang), minor_r(min_r), pitch(p) {};
+
+        vec3 integrate(vec3 (*integrand)(vec3 s, vec3 ds), scalar delta) const;
+
+        const char *name() const { return "Toroid"; };
+
+        friend std::istream operator>>(std::istream &is, LineSegment &ls);
+    };
+
+    std::istream operator>>(std::istream &is, LineSegment &ls);
+}
+#endif
diff --git a/include/fml/fml.h b/include/fml/fml.h
index 046db5a..bef1851 100644
--- a/include/fml/fml.h
+++ b/include/fml/fml.h
@@ -4,7 +4,16 @@
  * Copyright (C) 2019 Franklin Wei
  */
 
+#ifndef FML_H
+#define FML_H
+
 typedef float scalar;
 
+#include "curve.h"
+#include "manifold.h"
 #include "quat.h"
+#include "surface.h"
+#include "vec2.h"
 #include "vec3.h"
+
+#endif
diff --git a/include/fml/manifold.h b/include/fml/manifold.h
new file mode 100644
index 0000000..4f3664f
--- /dev/null
+++ b/include/fml/manifold.h
@@ -0,0 +1,20 @@
+#ifndef MANIFOLD_H
+#define MANIFOLD_H
+
+#include <cmath>
+#include <iostream>
+
+#include "fml.h"
+#include "vec3.h"
+
+namespace fml {
+    /* All manifolds inherit this class */
+    class Manifold {
+    public:
+        virtual vec3 integrate(vec3 (*integrand)(vec3 s, vec3 ds), scalar delta) const = 0;
+        virtual const char *name() const = 0;
+        virtual const int dimension() const = 0; // 0 = point, 1 = curve, 2 = surface, 3 = solid
+    };
+};
+
+#endif
diff --git a/include/fml/surface.h b/include/fml/surface.h
new file mode 100644
index 0000000..e42efe0
--- /dev/null
+++ b/include/fml/surface.h
@@ -0,0 +1,121 @@
+#ifndef SURFACE_H
+#define SURFACE_H
+
+#include "fml.h"
+#include "manifold.h"
+
+namespace fml {
+    /* All surfaces inherit this class */
+    /* The exact meaning of d is surface-dependent (see class
+     * definitions below); the limit of the calculated integral must
+     * approach its true value as d->0. d must be > 0. */
+    class Surface : public Manifold {
+    public:
+        const int dimension() const { return 2; }
+    };
+
+    class Plane : public Surface {
+    private:
+        /* The surface is specified by all points p = p0 + s v1 + t v2,
+         * such that 0 <= {s, t} < 1
+         *
+         * v1 and v2 must NOT be parallel.
+         *
+         * v1 and v2 should (but do not have to be) be normal (this is
+         * motivated primarily by usability; if the two vectors are indeed
+         * normal, then m1 and m2 have the nice geometric meaning of side
+         * length).
+         *
+         * d = ds = dt.
+         * dA will be in the direction of v1 x v2.
+         */
+        vec3 p0, v1, v2;
+
+    public:
+        Plane(vec3 p, vec3 _v1, vec3 _v2) : p0(p), v1(_v1), v2(_v2) {};
+
+        vec3 integrate(vec3 (*integrand)(vec3 s, vec3 dA), scalar d) const;
+        const char *name() const { return "Plane"; }
+    };
+
+/* Flat, circular disk (of varying extent) */
+    class Disk : public Surface {
+    private:
+        /* This represents a (possibly incomplete) circular disk in space,
+         * with a center, radius, and normal vector as specified.
+         *
+         * `angle' specifies the extent of the disk; angle=2pi for a full
+         * circle.
+         *
+         * We use:
+         * d = dr.
+         * d_theta = d / r.
+         *
+         * This makes it so that differential area elements on the edge of
+         * the disk are square.
+         *
+         * dA will be in the direction of normal.
+         */
+        vec3 center, radius, normal;
+        scalar angle;
+
+    public:
+        Disk(vec3 c, vec3 r, vec3 n, scalar a) : center(c), radius(r), normal(n), angle(a) {};
+
+        vec3 integrate(vec3 (*integrand)(vec3 s, vec3 dA), scalar d) const;
+        const char *name() const { return "Disk"; }
+    };
+
+/* Hollow, spherical shell */
+    class Sphere : public Surface {
+    private:
+        /* d = dtheta = dphi */
+
+        vec3 center;
+        scalar radius;
+
+    public:
+        Sphere(vec3 c, scalar r) : center(c), radius(r) {};
+
+        vec3 integrate(vec3 (*integrand)(vec3 s, vec3 dA), scalar d) const;
+        const char *name() const { return "Sphere"; }
+    };
+
+/* Cylinder without end caps */
+    class OpenCylinder : public Surface {
+    private:
+        /* Like this:
+         *
+         *         ___________________
+         *        / \                 \
+         *       origin ---------------> axis
+         *        \_/_________________/
+         *
+         */
+
+        vec3 origin, axis;
+        scalar radius;
+
+    public:
+        OpenCylinder(vec3 o, vec3 a, scalar r) : origin(o), axis(a), radius(r) {}
+
+        vec3 integrate(vec3 (*integrand)(vec3 s, vec3 dA), scalar d) const;
+        const char *name() const { return "OpenCylinder"; }
+    };
+
+/* Capped cylinder */
+    class ClosedCylinder : public OpenCylinder {
+    private:
+        Disk cap1, cap2;
+
+    public:
+        ClosedCylinder(vec3 o, vec3 a, scalar r) : OpenCylinder(o, a, r),
+                                                   cap1(o, vec3::any_unit_normal(a), -a.normalize(), 2*M_PI),
+                                                   cap2(o + a, vec3::any_unit_normal(a), a.normalize(), 2*M_PI) {}
+
+        vec3 integrate(vec3 (*integrand)(vec3 s, vec3 dA), scalar d) const;
+        const char *name() const { return "ClosedCylinder"; }
+    };
+}
+
+#endif