ODE physics
authorparasti <parasti@78b8d119-cf0a-0410-b17c-f493084dd1d7>
Mon, 11 Jan 2010 21:14:19 +0000 (21:14 +0000)
committerparasti <parasti@78b8d119-cf0a-0410-b17c-f493084dd1d7>
Mon, 11 Jan 2010 21:14:19 +0000 (21:14 +0000)
Disabled by default for now, enabled by compiling with "make
ENABLE_ODE=1".

Minimum ODE version is 0.11 due to use of kinematic bodies, but up to
the current version (which is 0.11.1) ODE also contains some sloppy
convex geom collision detection, so a patched version of ODE 0.11.1 is
available from

    http://github.com/parasti/ode/archives/master

Basically, no other ODE version will do.  Compile with "./configure
--enable-shared; make", install as usual.

Off the top of my head, things that don't work yet (and hopefully can
be made to):

 * Neverputt friction; Putt is pretty much broken.
 * Ball pendulums.
 * Displayed platform motion is out of sync with the simulation.
 * Bounces that land on the seams between lumps are not behaving well
   This is especially bad on curved surfaces, the ball's trajectory
   can be totally disrupted out of nowhere.  Cause might be solved in
   ODE code again.

git-svn-id: https://s.snth.net/svn/neverball/trunk@3119 78b8d119-cf0a-0410-b17c-f493084dd1d7

16 files changed:
Makefile
ball/game_server.c
putt/game.c
share/mapc.c
share/solid.c
share/solid.h
share/solid_all.c [new file with mode: 0644]
share/solid_all.h [new file with mode: 0644]
share/solid_cmd.c [new file with mode: 0644]
share/solid_cmd.h [new file with mode: 0644]
share/solid_gl.c
share/solid_phys.c [deleted file]
share/solid_phys.h [deleted file]
share/solid_sim.h [new file with mode: 0644]
share/solid_sim_ode.c [new file with mode: 0644]
share/solid_sim_sol.c [new file with mode: 0644]

index b45bcb0..d2d8640 100644 (file)
--- a/Makefile
+++ b/Makefile
@@ -87,6 +87,10 @@ ifeq ($(ENABLE_WII),1)
     ALL_CPPFLAGS += -DENABLE_WII=1
 endif
 
+ifeq ($(ENABLE_ODE),1)
+    ALL_CPPFLAGS += $(shell ode-config --cflags)
+endif
+
 ifdef DARWIN
     ALL_CPPFLAGS += -I/opt/local/include
 endif
@@ -111,6 +115,10 @@ ifeq ($(ENABLE_WII),1)
     TILT_LIBS := -lcwiimote -lbluetooth
 endif
 
+ifeq ($(ENABLE_ODE),1)
+    ODE_LIBS := $(shell ode-config --libs)
+endif
+
 OGL_LIBS := -lGL -lm
 
 ifdef MINGW
@@ -118,6 +126,10 @@ ifdef MINGW
         INTL_LIBS := -lintl
     endif
 
+    ifeq ($(ENABLE_ODE),1)
+        ODE_LIBS := -lode
+    endif
+
     TILT_LIBS :=
     OGL_LIBS  := -lopengl32 -lm
 endif
@@ -127,6 +139,10 @@ ifdef DARWIN
         INTL_LIBS := -lintl
     endif
 
+    ifeq ($(ENABLE_ODE),1)
+        ODE_LIBS := -lode
+    endif
+
     TILT_LIBS :=
     OGL_LIBS  := -framework OpenGL
 endif
@@ -138,7 +154,7 @@ ifdef DARWIN
 endif
 
 ALL_LIBS := $(SDL_LIBS) $(BASE_LIBS) $(TILT_LIBS) $(INTL_LIBS) -lSDL_ttf \
-    -lvorbisfile $(OGL_LIBS)
+    -lvorbisfile $(OGL_LIBS) $(ODE_LIBS)
 
 #------------------------------------------------------------------------------
 
@@ -180,7 +196,8 @@ BALL_OBJS := \
        share/image.o       \
        share/solid.o       \
        share/solid_gl.o    \
-       share/solid_phys.o  \
+       share/solid_cmd.o   \
+       share/solid_all.o   \
        share/part.o        \
        share/back.o        \
        share/geom.o        \
@@ -248,7 +265,8 @@ PUTT_OBJS := \
        share/image.o       \
        share/solid.o       \
        share/solid_gl.o    \
-       share/solid_phys.o  \
+       share/solid_cmd.o   \
+       share/solid_all.o   \
        share/part.o        \
        share/geom.o        \
        share/ball.o        \
@@ -280,6 +298,14 @@ PUTT_OBJS := \
        putt/st_conf.o      \
        putt/main.o
 
+ifeq ($(ENABLE_ODE),1)
+BALL_OBJS += share/solid_sim_ode.o
+PUTT_OBJS += share/solid_sim_ode.o
+else
+BALL_OBJS += share/solid_sim_sol.o
+PUTT_OBJS += share/solid_sim_sol.o
+endif
+
 ifdef MINGW
 BALL_OBJS += neverball.ico.o
 PUTT_OBJS += neverputt.ico.o
index b51fe60..12f97b3 100644 (file)
 
 #include "vec3.h"
 #include "item.h"
-#include "solid_phys.h"
 #include "config.h"
 #include "binary.h"
 #include "common.h"
 
+#include "solid_sim.h"
+#include "solid_all.h"
+#include "solid_cmd.h"
+
 #include "game_common.h"
 #include "game_server.h"
 #include "game_proxy.h"
@@ -543,6 +546,10 @@ int game_server_init(const char *file_name, int t, int e)
     got_orig = 0;
     grow = 0;
 
+    /* Initialize simulation. */
+
+    sol_init_sim(&file);
+
     sol_cmd_enq_func(game_proxy_enq);
 
     /* Queue client commands. */
@@ -563,6 +570,7 @@ void game_server_free(void)
 {
     if (server_state)
     {
+        sol_quit_sim();
         sol_free(&file);
         server_state = 0;
     }
index 92959e7..219ae48 100644 (file)
 #include "hud.h"
 #include "image.h"
 #include "audio.h"
-#include "solid_gl.h"
-#include "solid_phys.h"
 #include "config.h"
 #include "video.h"
 
+#include "solid_gl.h"
+#include "solid_sim.h"
+#include "solid_all.h"
+
 /*---------------------------------------------------------------------------*/
 
 static struct s_file file;
@@ -87,10 +89,12 @@ void game_init(const char *s)
 
     view_init();
     sol_load_gl(&file, s, config_get_d(CONFIG_SHADOW));
+    sol_init_sim(&file);
 }
 
 void game_free(void)
 {
+    sol_quit_sim();
     sol_free_gl(&file);
 }
 
index 41d2151..289f057 100644 (file)
@@ -1434,12 +1434,13 @@ static void clip_edge(struct s_file *fp,
 /*
  * Find all verts that lie on  the given side of the lump.  Sort these
  * verts to  have a counter-clockwise winding about  the plane normal.
- * Create geoms to tessellate the resulting convex polygon.
+ * Add the resulting convex polygon to the lump.
  */
-static void clip_geom(struct s_file *fp,
+static void clip_face(struct s_file *fp,
                       struct s_lump *lp, int si)
 {
-    int   m[256], t[256], d, i, j, n = 0;
+    int m[256], d, i, j, n = 0;
+
     float u[3];
     float v[3];
     float w[3];
@@ -1449,46 +1450,65 @@ static void clip_geom(struct s_file *fp,
     /* Find em. */
 
     for (i = 0; i < lp->vc; i++)
-    {
-        int vi = fp->iv[lp->v0 + i];
-
-        if (on_side(fp->vv[vi].p, sp))
-        {
-            m[n] = vi;
-            t[n] = inct(fp);
-
-            v_add(v, fp->vv[vi].p, plane_p[si]);
-
-            fp->tv[t[n]].u[0] = v_dot(v, plane_u[si]);
-            fp->tv[t[n]].u[1] = v_dot(v, plane_v[si]);
-
-            n++;
-        }
-    }
+        if (on_side(fp->vv[fp->iv[lp->v0 + i]].p, sp))
+            m[n++] = i;
 
     /* Sort em. */
 
     for (i = 1; i < n; i++)
         for (j = i + 1; j < n; j++)
         {
-            v_sub(u, fp->vv[m[i]].p, fp->vv[m[0]].p);
-            v_sub(v, fp->vv[m[j]].p, fp->vv[m[0]].p);
+            float *p0 = fp->vv[fp->iv[lp->v0 + m[0]]].p;
+            float *p1 = fp->vv[fp->iv[lp->v0 + m[i]]].p;
+            float *p2 = fp->vv[fp->iv[lp->v0 + m[j]]].p;
+
+            v_sub(u, p1, p0);
+            v_sub(v, p2, p0);
             v_crs(w, u, v);
 
-            if (v_dot(w, sp->n) < 0.f)
+            if (v_dot(w, sp->n) < 0.0f)
             {
-                d     = m[i];
-                m[i]  = m[j];
-                m[j]  =    d;
-
-                d     = t[i];
-                t[i]  = t[j];
-                t[j]  =    d;
+                d    = m[i];
+                m[i] = m[j];
+                m[j] = d;
             }
         }
 
     /* Index em. */
 
+    fp->iv[inci(fp)] = n;
+    lp->fc++;
+
+    for (i = 0; i < n; i++)
+    {
+        fp->iv[inci(fp)] = m[i];
+        lp->fc++;
+    }
+}
+
+/*
+ * Create geoms to tessellate the given convex polygon.
+ */
+static void clip_geom(struct s_file *fp,
+                      struct s_lump *lp, int si, int fi)
+{
+    int   t[256], i, n;
+    float v[3];
+
+    n = fp->iv[fi++];
+
+    for (i = 0; i < n; i++)
+    {
+        int vi = fp->iv[lp->v0 + fp->iv[fi + i]];
+
+        t[i] = inct(fp);
+
+        v_add(v, fp->vv[vi].p, plane_p[si]);
+
+        fp->tv[t[i]].u[0] = v_dot(v, plane_u[si]);
+        fp->tv[t[i]].u[1] = v_dot(v, plane_v[si]);
+    }
+
     for (i = 0; i < n - 2; i++)
     {
         fp->gv[fp->gc].mi = plane_m[si];
@@ -1501,9 +1521,9 @@ static void clip_geom(struct s_file *fp,
         fp->gv[fp->gc].sj = si;
         fp->gv[fp->gc].sk = si;
 
-        fp->gv[fp->gc].vi = m[0];
-        fp->gv[fp->gc].vj = m[i + 1];
-        fp->gv[fp->gc].vk = m[i + 2];
+        fp->gv[fp->gc].vi = fp->iv[lp->v0 + fp->iv[fi]];
+        fp->gv[fp->gc].vj = fp->iv[lp->v0 + fp->iv[fi + i + 1]];
+        fp->gv[fp->gc].vk = fp->iv[lp->v0 + fp->iv[fi + i + 2]];
 
         fp->iv[fp->ic] = fp->gc;
         inci(fp);
@@ -1519,7 +1539,7 @@ static void clip_geom(struct s_file *fp,
  */
 static void clip_lump(struct s_file *fp, struct s_lump *lp)
 {
-    int i, j, k;
+    int i, j, k, fi;
 
     lp->v0 = fp->ic;
     lp->vc = 0;
@@ -1541,13 +1561,24 @@ static void clip_lump(struct s_file *fp, struct s_lump *lp)
                       fp->iv[lp->s0 + i],
                       fp->iv[lp->s0 + j]);
 
+    lp->f0 = fp->ic;
+    lp->fc = 0;
+
+    for (i = 0; i < lp->sc; i++)
+        clip_face(fp, lp, fp->iv[lp->s0 + i]);
+
     lp->g0 = fp->ic;
     lp->gc = 0;
 
+    fi = lp->f0;
+
     for (i = 0; i < lp->sc; i++)
+    {
         if (fp->mv[plane_m[fp->iv[lp->s0 + i]]].d[3] > 0.0f)
-            clip_geom(fp, lp,
-                      fp->iv[lp->s0 + i]);
+            clip_geom(fp, lp, fp->iv[lp->s0 + i], fi);
+
+        fi += fp->iv[fi] + 1;
+    }
 
     for (i = 0; i < lp->sc; i++)
         if (plane_f[fp->iv[lp->s0 + i]])
index 9a8e170..d080e33 100644 (file)
@@ -22,7 +22,7 @@
 #include "fs.h"
 
 #define MAGIC       0x4c4f53af
-#define SOL_VERSION 6
+#define SOL_VERSION 7
 
 /*---------------------------------------------------------------------------*/
 
@@ -85,6 +85,8 @@ static void sol_load_lump(fs_file fin, struct s_lump *lp)
     get_index(fin, &lp->gc);
     get_index(fin, &lp->s0);
     get_index(fin, &lp->sc);
+    get_index(fin, &lp->f0);
+    get_index(fin, &lp->fc);
 }
 
 static void sol_load_node(fs_file fin, struct s_node *np)
@@ -433,6 +435,8 @@ static void sol_stor_lump(fs_file fout, struct s_lump *lp)
     put_index(fout, &lp->gc);
     put_index(fout, &lp->s0);
     put_index(fout, &lp->sc);
+    put_index(fout, &lp->f0);
+    put_index(fout, &lp->fc);
 }
 
 static void sol_stor_node(fs_file fout, struct s_node *np)
index 7f5d643..0d682c6 100644 (file)
@@ -160,6 +160,7 @@ struct s_lump
     int e0, ec;
     int g0, gc;
     int s0, sc;
+    int f0, fc;                                /* polygon data               */
 };
 
 struct s_node
diff --git a/share/solid_all.c b/share/solid_all.c
new file mode 100644 (file)
index 0000000..98af0d5
--- /dev/null
@@ -0,0 +1,467 @@
+
+/* Random code used in more than one place. */
+
+#include "solid_all.h"
+#include "solid_cmd.h"
+#include "solid.h"
+#include "common.h"
+#include "vec3.h"
+#include "geom.h"
+
+/*---------------------------------------------------------------------------*/
+
+static float erp(float t)
+{
+    return 3.0f * t * t - 2.0f * t * t * t;
+}
+
+#if UNUSED
+static float derp(float t)
+{
+    return 6.0f * t     - 6.0f * t * t;
+}
+#endif
+
+void sol_body_p(float p[3], const struct s_file *fp, int pi, float t)
+{
+    float v[3];
+
+    if (pi >= 0)
+    {
+        const struct s_path *pp = fp->pv + pi;
+        const struct s_path *pq = fp->pv + pp->pi;
+
+        float s = MIN(t / pp->t, 1.0f);
+
+        v_sub(v, pq->p, pp->p);
+        v_mad(p, pp->p, v, pp->s ? erp(s) : s);
+    }
+    else
+    {
+        p[0] = 0.0f;
+        p[1] = 0.0f;
+        p[2] = 0.0f;
+    }
+}
+
+void sol_body_v(float v[3],
+                const struct s_file *fp,
+                int pi, float t, float dt)
+{
+    if (pi >= 0 && fp->pv[pi].f)
+    {
+        float p[3], q[3];
+
+        sol_body_p(p, fp, pi, t);
+        sol_body_p(q, fp, pi, t + dt);
+
+        v_sub(v, q, p);
+
+        v[0] /= dt;
+        v[1] /= dt;
+        v[2] /= dt;
+    }
+    else
+    {
+        v[0] = 0.0f;
+        v[1] = 0.0f;
+        v[2] = 0.0f;
+    }
+}
+
+/*---------------------------------------------------------------------------*/
+
+/*
+ * Integrate the rotation of the given basis E under angular velocity W
+ * through time DT.
+ */
+void sol_rotate(float e[3][3], const float w[3], float dt)
+{
+    if (v_len(w) > 0.0f)
+    {
+        float a[3], M[16], f[3][3];
+
+        /* Compute the rotation matrix. */
+
+        v_nrm(a, w);
+        m_rot(M, a, v_len(w) * dt);
+
+        /* Apply it to the basis. */
+
+        m_vxfm(f[0], M, e[0]);
+        m_vxfm(f[1], M, e[1]);
+        m_vxfm(f[2], M, e[2]);
+
+        /* Re-orthonormalize the basis. */
+
+        v_crs(e[2], f[0], f[1]);
+        v_crs(e[1], f[2], f[0]);
+        v_crs(e[0], f[1], f[2]);
+
+        v_nrm(e[0], e[0]);
+        v_nrm(e[1], e[1]);
+        v_nrm(e[2], e[2]);
+    }
+}
+
+/*---------------------------------------------------------------------------*/
+
+/*
+ * Compute the new angular velocity and orientation of a ball pendulum.
+ * A gives the accelleration of the ball.  G gives the gravity vector.
+ */
+void sol_pendulum(struct s_ball *up,
+                  const float a[3],
+                  const float g[3], float dt)
+{
+    float v[3], A[3], F[3], r[3], Y[3], T[3] = { 0.0f, 0.0f, 0.0f };
+
+    const float m  = 5.000f;
+    const float ka = 0.500f;
+    const float kd = 0.995f;
+
+    /* Find the total force over DT. */
+
+    v_scl(A, a,     ka);
+    v_mad(A, A, g, -dt);
+
+    /* Find the force. */
+
+    v_scl(F, A, m / dt);
+
+    /* Find the position of the pendulum. */
+
+    v_scl(r, up->E[1], -up->r);
+
+    /* Find the torque on the pendulum. */
+
+    if (fabsf(v_dot(r, F)) > 0.0f)
+        v_crs(T, F, r);
+
+    /* Apply the torque and dampen the angular velocity. */
+
+    v_mad(up->W, up->W, T, dt);
+    v_scl(up->W, up->W,    kd);
+
+    /* Apply the angular velocity to the pendulum basis. */
+
+    sol_rotate(up->E, up->W, dt);
+
+    /* Apply a torque turning the pendulum toward the ball velocity. */
+
+    v_mad(v, up->v, up->E[1], v_dot(up->v, up->E[1]));
+    v_crs(Y, v, up->E[2]);
+    v_scl(Y, up->E[1], 2 * v_dot(Y, up->E[1]));
+
+    sol_rotate(up->E, Y, dt);
+}
+
+/*---------------------------------------------------------------------------*/
+
+/*
+ * Compute the states of all switches after DT seconds have passed.
+ */
+void sol_swch_step(struct s_file *fp, float dt)
+{
+    int xi;
+
+    union cmd cmd;
+
+    for (xi = 0; xi < fp->xc; xi++)
+    {
+        struct s_swch *xp = fp->xv + xi;
+
+        volatile float t = xp->t;
+
+        if (t < xp->t0)
+        {
+            xp->t = (t += dt);
+
+            if (t >= xp->t0)
+            {
+                int pi = xp->pi;
+                int pj = xp->pi;
+
+                do  /* Tortoise and hare cycle traverser. */
+                {
+                    fp->pv[pi].f = xp->f0;
+                    fp->pv[pj].f = xp->f0;
+
+                    cmd.type        = CMD_PATH_FLAG;
+                    cmd.pathflag.pi = pi;
+                    cmd.pathflag.f  = fp->pv[pi].f;
+                    sol_cmd_enq(&cmd);
+
+                    pi = fp->pv[pi].pi;
+                    pj = fp->pv[pj].pi;
+                    pj = fp->pv[pj].pi;
+                }
+                while (pi != pj);
+
+                xp->f = xp->f0;
+
+                cmd.type          = CMD_SWCH_TOGGLE;
+                cmd.swchtoggle.xi = xi;
+                sol_cmd_enq(&cmd);
+            }
+        }
+    }
+}
+
+/*
+ * Compute the positions of all bodies after DT seconds have passed.
+ */
+void sol_body_step(struct s_file *fp, float dt)
+{
+    int i;
+
+    union cmd cmd;
+
+    for (i = 0; i < fp->bc; i++)
+    {
+        struct s_body *bp = fp->bv + i;
+        struct s_path *pp = fp->pv + bp->pi;
+
+        volatile float t = bp->t;
+
+        if (bp->pi >= 0 && pp->f)
+        {
+            bp->t = (t += dt);
+
+            if (t >= pp->t)
+            {
+                bp->t  = 0;
+                bp->pi = pp->pi;
+
+                cmd.type        = CMD_BODY_TIME;
+                cmd.bodytime.bi = i;
+                cmd.bodytime.t  = bp->t;
+                sol_cmd_enq(&cmd);
+
+                cmd.type        = CMD_BODY_PATH;
+                cmd.bodypath.bi = i;
+                cmd.bodypath.pi = bp->pi;
+                sol_cmd_enq(&cmd);
+            }
+        }
+    }
+}
+
+/*
+ * Compute the positions of all balls after DT seconds have passed.
+ */
+void sol_ball_step(struct s_file *fp, float dt)
+{
+    int i;
+
+    for (i = 0; i < fp->uc; i++)
+    {
+        struct s_ball *up = fp->uv + i;
+
+        v_mad(up->p, up->p, up->v, dt);
+
+        sol_rotate(up->e, up->w, dt);
+    }
+}
+
+/*---------------------------------------------------------------------------*/
+
+int sol_item_test(struct s_file *fp, float *p, float item_r)
+{
+    const float *ball_p = fp->uv->p;
+    const float  ball_r = fp->uv->r;
+
+    int hi;
+
+    for (hi = 0; hi < fp->hc; hi++)
+    {
+        float r[3];
+
+        r[0] = ball_p[0] - fp->hv[hi].p[0];
+        r[1] = ball_p[1] - fp->hv[hi].p[1];
+        r[2] = ball_p[2] - fp->hv[hi].p[2];
+
+        if (fp->hv[hi].t != ITEM_NONE && v_len(r) < ball_r + item_r)
+        {
+            p[0] = fp->hv[hi].p[0];
+            p[1] = fp->hv[hi].p[1];
+            p[2] = fp->hv[hi].p[2];
+
+            return hi;
+        }
+    }
+    return -1;
+}
+
+struct s_goal *sol_goal_test(struct s_file *fp, float *p, int ui)
+{
+    const float *ball_p = fp->uv[ui].p;
+    const float  ball_r = fp->uv[ui].r;
+    int zi;
+
+    for (zi = 0; zi < fp->zc; zi++)
+    {
+        float r[3];
+
+        r[0] = ball_p[0] - fp->zv[zi].p[0];
+        r[1] = ball_p[2] - fp->zv[zi].p[2];
+        r[2] = 0;
+
+        if (v_len(r) < fp->zv[zi].r - ball_r &&
+            ball_p[1] > fp->zv[zi].p[1] &&
+            ball_p[1] < fp->zv[zi].p[1] + GOAL_HEIGHT / 2)
+        {
+            p[0] = fp->zv[zi].p[0];
+            p[1] = fp->zv[zi].p[1];
+            p[2] = fp->zv[zi].p[2];
+
+            return &fp->zv[zi];
+        }
+    }
+    return NULL;
+}
+
+/*
+ * Test if the  ball UI is inside a  jump. Return 1 if yes  and fill P
+ * with the destination position, return 0 if not, and return 2 if the
+ * ball is on the border of a jump.
+ */
+int sol_jump_test(struct s_file *fp, float *p, int ui)
+{
+    const float *ball_p = fp->uv[ui].p;
+    const float  ball_r = fp->uv[ui].r;
+    int ji;
+    float l;
+    int res = 0;
+
+    for (ji = 0; ji < fp->jc; ji++)
+    {
+        float r[3];
+
+        r[0] = ball_p[0] - fp->jv[ji].p[0];
+        r[1] = ball_p[2] - fp->jv[ji].p[2];
+        r[2] = 0;
+
+        l = v_len(r) - fp->jv[ji].r;
+        if (l < 0 &&
+            ball_p[1] > fp->jv[ji].p[1] &&
+            ball_p[1] < fp->jv[ji].p[1] + JUMP_HEIGHT / 2)
+        {
+            if (l < - ball_r )
+            {
+                p[0] = fp->jv[ji].q[0] + (ball_p[0] - fp->jv[ji].p[0]);
+                p[1] = fp->jv[ji].q[1] + (ball_p[1] - fp->jv[ji].p[1]);
+                p[2] = fp->jv[ji].q[2] + (ball_p[2] - fp->jv[ji].p[2]);
+
+                return 1;
+            }
+            else
+                res = 2;
+        }
+    }
+    return res;
+}
+
+/*
+ * Test and process the event the ball UI enters a switch. Return 1 if
+ * a visible  switch is  activated, return 0  otherwise (no  switch is
+ * activated or only invisible switches).
+ */
+int sol_swch_test(struct s_file *fp, int ui)
+{
+    const float *ball_p = fp->uv[ui].p;
+    const float  ball_r = fp->uv[ui].r;
+    int xi;
+    int res = 0;
+
+    union cmd cmd;
+
+    for (xi = 0; xi < fp->xc; xi++)
+    {
+        struct s_swch *xp = fp->xv + xi;
+
+        /* FIXME enter/exit events don't work for timed switches */
+
+        if (xp->t0 == 0 || xp->f == xp->f0)
+        {
+            float l;
+            float r[3];
+
+            r[0] = ball_p[0] - xp->p[0];
+            r[1] = ball_p[2] - xp->p[2];
+            r[2] = 0;
+
+            l = v_len(r) - xp->r;
+
+            if (l < ball_r &&
+                ball_p[1] > xp->p[1] &&
+                ball_p[1] < xp->p[1] + SWCH_HEIGHT / 2)
+            {
+                if (!xp->e && l < - ball_r)
+                {
+                    int pi = xp->pi;
+                    int pj = xp->pi;
+
+                    /* The ball enters. */
+
+                    if (xp->t0 == 0)
+                    {
+                        xp->e = 1;
+
+                        cmd.type         = CMD_SWCH_ENTER;
+                        cmd.swchenter.xi = xi;
+                        sol_cmd_enq(&cmd);
+                    }
+
+                    /* Toggle the state, update the path. */
+
+                    xp->f = xp->f ? 0 : 1;
+
+                    cmd.type          = CMD_SWCH_TOGGLE;
+                    cmd.swchtoggle.xi = xi;
+                    sol_cmd_enq(&cmd);
+
+                    do  /* Tortoise and hare cycle traverser. */
+                    {
+                        fp->pv[pi].f = xp->f;
+                        fp->pv[pj].f = xp->f;
+
+                        cmd.type        = CMD_PATH_FLAG;
+                        cmd.pathflag.pi = pi;
+                        cmd.pathflag.f  = fp->pv[pi].f;
+                        sol_cmd_enq(&cmd);
+
+                        pi = fp->pv[pi].pi;
+                        pj = fp->pv[pj].pi;
+                        pj = fp->pv[pj].pi;
+                    }
+                    while (pi != pj);
+
+                    /* It toggled to non-default state, start the timer. */
+
+                    if (xp->f != xp->f0)
+                        xp->t = 0.0f;
+
+                    /* If visible, set the result. */
+
+                    if (!xp->i)
+                        res = 1;
+                }
+            }
+
+            /* The ball exits. */
+
+            else if (xp->e)
+            {
+                xp->e = 0;
+
+                cmd.type        = CMD_SWCH_EXIT;
+                cmd.swchexit.xi = xi;
+                sol_cmd_enq(&cmd);
+            }
+        }
+    }
+    return res;
+}
+
+/*---------------------------------------------------------------------------*/
diff --git a/share/solid_all.h b/share/solid_all.h
new file mode 100644 (file)
index 0000000..17c5939
--- /dev/null
@@ -0,0 +1,24 @@
+#ifndef SOLID_ALL_H
+#define SOLID_ALL_H
+
+#include "solid.h"
+
+void sol_body_p(float p[3], const struct s_file *fp, int pi, float t);
+void sol_body_v(float v[3], const struct s_file *fp, int pi, float t, float dt);
+
+void sol_rotate(float e[3][3], const float w[3], float dt);
+
+void sol_pendulum(struct s_ball *up,
+                  const float a[3],
+                  const float g[3], float dt);
+
+void sol_swch_step(struct s_file *fp, float dt);
+void sol_body_step(struct s_file *fp, float dt);
+void sol_ball_step(struct s_file *fp, float dt);
+
+int            sol_item_test(struct s_file *fp, float *p, float item_r);
+struct s_goal *sol_goal_test(struct s_file *fp, float *p, int ui);
+int            sol_jump_test(struct s_file *fp, float *p, int ui);
+int            sol_swch_test(struct s_file *fp, int ui);
+
+#endif
diff --git a/share/solid_cmd.c b/share/solid_cmd.c
new file mode 100644 (file)
index 0000000..c3f3316
--- /dev/null
@@ -0,0 +1,100 @@
+#include <stdlib.h>
+
+#include "solid_cmd.h"
+#include "cmd.h"
+#include "list.h"
+
+/*---------------------------------------------------------------------------*/
+
+/*
+ * sol_step generates a few commands that supersede previous commands
+ * of the same type generated by the same sol_step invocation.  The
+ * code below provides the means to accumulate commands for later
+ * addition to the command queue, while making sure that commands that
+ * have been superseded during the accumulation are discarded.
+ */
+
+int sol_cmd_defer;
+
+static void (*cmd_enq_fn)(const union cmd *);
+static List deferred_cmds;
+
+void sol_cmd_enq_deferred(void)
+{
+    List l, r;
+
+    /* Reverse the list to preserve original command order. */
+
+    for (r = NULL, l = deferred_cmds;
+         l;
+         r = list_cons(l->data, r), l = list_rest(l));
+
+    /* Enqueue commands. */
+
+    for (; r; r = list_rest(r))
+    {
+        if (cmd_enq_fn)
+            cmd_enq_fn(r->data);
+
+        free(r->data);
+    }
+
+    deferred_cmds = NULL;
+}
+
+void sol_cmd_enq(const union cmd *new)
+{
+    if (sol_cmd_defer)
+    {
+        union cmd *copy;
+        List l, p;
+
+        for (p = NULL, l = deferred_cmds; l; p = l, l = l->next)
+        {
+            union cmd *cur = l->data;
+
+            /* Remove element made obsolete by the new command. */
+
+            if (new->type == cur->type &&
+                ((new->type == CMD_BODY_TIME &&
+                  new->bodytime.bi == cur->bodytime.bi) ||
+                 (new->type == CMD_BODY_PATH &&
+                  new->bodypath.bi == cur->bodypath.bi)))
+            {
+                free(cur);
+
+                if (p)
+                    p->next       = list_rest(l);
+                else
+                    deferred_cmds = list_rest(l);
+
+                /*
+                 * The operation above made the list pointer useless
+                 * for the variable update part of the loop, and it
+                 * seems a bit involved to recover from that in a
+                 * proper fashion.  Fortunately, this very block
+                 * ensures that there's only one element to remove, so
+                 * no more iterations are needed.
+                 */
+
+                l = NULL;
+                break;
+            }
+        }
+
+        if ((copy = malloc(sizeof (*copy))))
+        {
+            *copy = *new;
+            deferred_cmds = list_cons(copy, deferred_cmds);
+        }
+    }
+    else if (cmd_enq_fn)
+        cmd_enq_fn(new);
+}
+
+void sol_cmd_enq_func(void (*enq_fn) (const union cmd *))
+{
+    cmd_enq_fn = enq_fn;
+}
+
+/*---------------------------------------------------------------------------*/
diff --git a/share/solid_cmd.h b/share/solid_cmd.h
new file mode 100644 (file)
index 0000000..1601d29
--- /dev/null
@@ -0,0 +1,12 @@
+#ifndef SOLID_CMD_H
+#define SOLID_CMD_H
+
+#include "cmd.h"
+
+extern int sol_cmd_defer;
+
+void sol_cmd_enq(const union cmd *);
+void sol_cmd_enq_deferred(void);
+void sol_cmd_enq_func(void (*enq_fn) (const union cmd *));
+
+#endif
index cd0282c..72214a4 100644 (file)
@@ -25,7 +25,7 @@
 #include "image.h"
 #include "base_image.h"
 #include "solid_gl.h"
-#include "solid_phys.h"
+#include "solid_all.h"
 #include "base_config.h"
 #include "lang.h"
 
diff --git a/share/solid_phys.c b/share/solid_phys.c
deleted file mode 100644 (file)
index 0f38df5..0000000
+++ /dev/null
@@ -1,1167 +0,0 @@
-/*
- * Copyright (C) 2003 Robert Kooima
- *
- * NEVERBALL is  free software; you can redistribute  it and/or modify
- * it under the  terms of the GNU General  Public License as published
- * by the Free  Software Foundation; either version 2  of the License,
- * or (at your option) any later version.
- *
- * This program is distributed in the hope that it will be useful, but
- * WITHOUT  ANY  WARRANTY;  without   even  the  implied  warranty  of
- * MERCHANTABILITY or  FITNESS FOR A PARTICULAR PURPOSE.   See the GNU
- * General Public License for more details.
- */
-
-#include <math.h>
-
-#include "vec3.h"
-#include "geom.h" /* Only for height constants! */
-#include "solid.h"
-#include "solid_phys.h"
-#include "common.h"
-
-#define LARGE 1.0e+5f
-#define SMALL 1.0e-3f
-
-/*---------------------------------------------------------------------------*/
-
-#include "cmd.h"
-#include "list.h"
-
-static union cmd cmd;
-
-static void (*cmd_enq_fn)(const union cmd *);
-
-/*
- * sol_step generates a few commands that supersede previous commands
- * of the same type generated by the same sol_step invocation.  The
- * code below provides the means to accumulate commands for later
- * addition to the command queue, while making sure that commands that
- * have been superseded during the accumulation are discarded.
- */
-
-static int  defer_cmds;
-static List deferred_cmds;
-
-static void sol_cmd_enq_deferred(void)
-{
-    List l, r;
-
-    /* Reverse the list to preserve original command order. */
-
-    for (r = NULL, l = deferred_cmds;
-         l;
-         r = list_cons(l->data, r), l = list_rest(l));
-
-    /* Enqueue commands. */
-
-    for (; r; r = list_rest(r))
-    {
-        if (cmd_enq_fn)
-            cmd_enq_fn(r->data);
-
-        free(r->data);
-    }
-
-    deferred_cmds = NULL;
-}
-
-static void sol_cmd_enq(const union cmd *new)
-{
-    if (defer_cmds)
-    {
-        union cmd *copy;
-        List l, p;
-
-        for (p = NULL, l = deferred_cmds; l; p = l, l = l->next)
-        {
-            union cmd *cur = l->data;
-
-            /* Remove element made obsolete by the new command. */
-
-            if (new->type == cur->type &&
-                ((new->type == CMD_BODY_TIME &&
-                  new->bodytime.bi == cur->bodytime.bi) ||
-                 (new->type == CMD_BODY_PATH &&
-                  new->bodypath.bi == cur->bodypath.bi)))
-            {
-                free(cur);
-
-                if (p)
-                    p->next       = list_rest(l);
-                else
-                    deferred_cmds = list_rest(l);
-
-                /*
-                 * The operation above made the list pointer useless
-                 * for the variable update part of the loop, and it
-                 * seems a bit involved to recover from that in a
-                 * proper fashion.  Fortunately, this very block
-                 * ensures that there's only one element to remove, so
-                 * no more iterations are needed.
-                 */
-
-                l = NULL;
-                break;
-            }
-        }
-
-        if ((copy = malloc(sizeof (*copy))))
-        {
-            *copy = *new;
-            deferred_cmds = list_cons(copy, deferred_cmds);
-        }
-    }
-    else if (cmd_enq_fn)
-        cmd_enq_fn(new);
-}
-
-void sol_cmd_enq_func(void (*enq_fn) (const union cmd *))
-{
-    cmd_enq_fn = enq_fn;
-}
-
-/*---------------------------------------------------------------------------*/
-
-static float erp(float t)
-{
-    return 3.0f * t * t - 2.0f * t * t * t;
-}
-
-#if UNUSED
-static float derp(float t)
-{
-    return 6.0f * t     - 6.0f * t * t;
-}
-#endif
-
-static void sol_body_v(float v[3],
-                       const struct s_file *fp,
-                       int pi, float t, float dt)
-{
-    if (pi >= 0 && fp->pv[pi].f)
-    {
-        float p[3], q[3];
-
-        sol_body_p(p, fp, pi, t);
-        sol_body_p(q, fp, pi, t + dt);
-
-        v_sub(v, q, p);
-
-        v[0] /= dt;
-        v[1] /= dt;
-        v[2] /= dt;
-    }
-    else
-    {
-        v[0] = 0.0f;
-        v[1] = 0.0f;
-        v[2] = 0.0f;
-    }
-}
-
-void sol_body_p(float p[3], const struct s_file *fp, int pi, float t)
-{
-    float v[3];
-
-    if (pi >= 0)
-    {
-        const struct s_path *pp = fp->pv + pi;
-        const struct s_path *pq = fp->pv + pp->pi;
-
-        float s = MIN(t / pp->t, 1.0f);
-
-        v_sub(v, pq->p, pp->p);
-        v_mad(p, pp->p, v, pp->s ? erp(s) : s);
-    }
-    else
-    {
-        p[0] = 0.0f;
-        p[1] = 0.0f;
-        p[2] = 0.0f;
-    }
-}
-
-/*---------------------------------------------------------------------------*/
-/* Solves (p + v * t) . (p + v * t) == r * r for smallest t.                 */
-
-/*
- * Given vectors A = P, B = V * t, C = A + B, |C| = r, solve for
- * smallest t.
- *
- * Some useful dot product properties:
- *
- * 1) A . A = |A| * |A|
- * 2) A . (B + C) = A . B + A . C
- * 3) (r * A) . B = r * (A . B)
- *
- * Deriving a quadratic equation:
- *
- * C . C = r * r                                     (1)
- * (A + B) . (A + B) = r * r
- * A . (A + B) + B . (A + B) = r * r                 (2)
- * A . A + A . B + B . A + B . B = r * r             (2)
- * A . A + 2 * (A . B) + B . B = r * r
- * P . P + 2 * (P . V * t) + (V * t . V * t) = r * r
- * P . P + 2 * (P . V) * t + (V . V) * t * t = r * r (3)
- * (V . V) * t * t + 2 * (P . V) * t + P . P - r * r = 0
- *
- * This equation is solved using the quadratic formula.
- */
-
-static float v_sol(const float p[3], const float v[3], float r)
-{
-    float a = v_dot(v, v);
-    float b = v_dot(v, p) * 2.0f;
-    float c = v_dot(p, p) - r * r;
-    float d = b * b - 4.0f * a * c;
-
-/* HACK: This seems to cause failures to detect low-velocity collision
-         Yet, the potential division by zero below seems fine.
-    if (fabsf(a) < SMALL) return LARGE;
-*/
-
-    if      (d < 0.0f) return LARGE;
-    else if (d > 0.0f)
-    {
-        float t0 = 0.5f * (-b - fsqrtf(d)) / a;
-        float t1 = 0.5f * (-b + fsqrtf(d)) / a;
-        float t  = (t0 < t1) ? t0 : t1;
-
-        return (t < 0.0f) ? LARGE : t;
-    }
-    else return -b * 0.5f / a;
-}
-
-/*---------------------------------------------------------------------------*/
-
-/*
- * Compute the  earliest time  and position of  the intersection  of a
- * sphere and a vertex.
- *
- * The sphere has radius R and moves along vector V from point P.  The
- * vertex moves  along vector  W from point  Q in a  coordinate system
- * based at O.
- */
-static float v_vert(float Q[3],
-                    const float o[3],
-                    const float q[3],
-                    const float w[3],
-                    const float p[3],
-                    const float v[3], float r)
-{
-    float O[3], P[3], V[3];
-    float t = LARGE;
-
-    v_add(O, o, q);
-    v_sub(P, p, O);
-    v_sub(V, v, w);
-
-    if (v_dot(P, V) < 0.0f)
-    {
-        t = v_sol(P, V, r);
-
-        if (t < LARGE)
-            v_mad(Q, O, w, t);
-    }
-    return t;
-}
-
-/*
- * Compute the  earliest time  and position of  the intersection  of a
- * sphere and an edge.
- *
- * The sphere has radius R and moves along vector V from point P.  The
- * edge moves along vector W from point Q in a coordinate system based
- * at O.  The edge extends along the length of vector U.
- */
-static float v_edge(float Q[3],
-                    const float o[3],
-                    const float q[3],
-                    const float u[3],
-                    const float w[3],
-                    const float p[3],
-                    const float v[3], float r)
-{
-    float d[3], e[3];
-    float P[3], V[3];
-    float du, eu, uu, s, t;
-
-    v_sub(d, p, o);
-    v_sub(d, d, q);
-    v_sub(e, v, w);
-
-    /*
-     * Think projections.  Vectors D, extending from the edge vertex Q
-     * to the sphere,  and E, the relative velocity  of sphere wrt the
-     * edge, are  made orthogonal to  the edge vector U.   Division of
-     * the  dot products  is required  to obtain  the  true projection
-     * ratios since U does not have unit length.
-     */
-
-    du = v_dot(d, u);
-    eu = v_dot(e, u);
-    uu = v_dot(u, u);
-
-    v_mad(P, d, u, -du / uu);
-    v_mad(V, e, u, -eu / uu);
-
-    t = v_sol(P, V, r);
-    s = (du + eu * t) / uu; /* Projection of D + E * t on U. */
-
-    if (0.0f <= t && t < LARGE && 0.0f < s && s < 1.0f)
-    {
-        v_mad(d, o, w, t);
-        v_mad(e, q, u, s);
-        v_add(Q, e, d);
-    }
-    else
-        t = LARGE;
-
-    return t;
-}
-
-/*
- * Compute  the earliest  time and  position of  the intersection  of a
- * sphere and a plane.
- *
- * The sphere has radius R and moves along vector V from point P.  The
- * plane  moves  along  vector  W.   The  plane has  normal  N  and  is
- * positioned at distance D from the origin O along that normal.
- */
-static float v_side(float Q[3],
-                    const float o[3],
-                    const float w[3],
-                    const float n[3], float d,
-                    const float p[3],
-                    const float v[3], float r)
-{
-    float vn = v_dot(v, n);
-    float wn = v_dot(w, n);
-    float t  = LARGE;
-
-    if (vn - wn <= 0.0f)
-    {
-        float on = v_dot(o, n);
-        float pn = v_dot(p, n);
-
-        float u = (r + d + on - pn) / (vn - wn);
-        float a = (    d + on - pn) / (vn - wn);
-
-        if (0.0f <= u)
-        {
-            t = u;
-
-            v_mad(Q, p, v, +t);
-            v_mad(Q, Q, n, -r);
-        }
-        else if (0.0f <= a)
-        {
-            t = 0;
-
-            v_mad(Q, p, v, +t);
-            v_mad(Q, Q, n, -r);
-        }
-    }
-    return t;
-}
-
-/*---------------------------------------------------------------------------*/
-
-/*
- * Integrate the rotation of the given basis E under angular velocity W
- * through time DT.
- */
-static void sol_rotate(float e[3][3], const float w[3], float dt)
-{
-    if (v_len(w) > 0.0f)
-    {
-        float a[3], M[16], f[3][3];
-
-        /* Compute the rotation matrix. */
-
-        v_nrm(a, w);
-        m_rot(M, a, v_len(w) * dt);
-
-        /* Apply it to the basis. */
-
-        m_vxfm(f[0], M, e[0]);
-        m_vxfm(f[1], M, e[1]);
-        m_vxfm(f[2], M, e[2]);
-
-        /* Re-orthonormalize the basis. */
-
-        v_crs(e[2], f[0], f[1]);
-        v_crs(e[1], f[2], f[0]);
-        v_crs(e[0], f[1], f[2]);
-
-        v_nrm(e[0], e[0]);
-        v_nrm(e[1], e[1]);
-        v_nrm(e[2], e[2]);
-    }
-}
-
-/*
- * Compute the new  linear and angular velocities of  a bouncing ball.
- * Q  gives the  position  of the  point  of impact  and  W gives  the
- * velocity of the object being impacted.
- */
-static float sol_bounce(struct s_ball *up,
-                        const float q[3],
-                        const float w[3], float dt)
-{
-    float n[3], r[3], d[3], vn, wn;
-    float *p = up->p;
-    float *v = up->v;
-
-    /* Find the normal of the impact. */
-
-    v_sub(r, p, q);
-    v_sub(d, v, w);
-    v_nrm(n, r);
-
-    /* Find the new angular velocity. */
-
-    v_crs(up->w, d, r);
-    v_scl(up->w, up->w, -1.0f / (up->r * up->r));
-
-    /* Find the new linear velocity. */
-
-    vn = v_dot(v, n);
-    wn = v_dot(w, n);
-
-    v_mad(v, v, n, 1.7 * (wn - vn));
-
-    v_mad(p, q, n, up->r);
-
-    /* Return the "energy" of the impact, to determine the sound amplitude. */
-
-    return fabsf(v_dot(n, d));
-}
-
-/*
- * Compute the new angular velocity and orientation of a ball pendulum.
- * A gives the accelleration of the ball.  G gives the gravity vector.
- */
-static void sol_pendulum(struct s_ball *up,
-                         const float a[3],
-                         const float g[3], float dt)
-{
-    float v[3], A[3], F[3], r[3], Y[3], T[3] = { 0.0f, 0.0f, 0.0f };
-
-    const float m  = 5.000f;
-    const float ka = 0.500f;
-    const float kd = 0.995f;
-
-    /* Find the total force over DT. */
-
-    v_scl(A, a,     ka);
-    v_mad(A, A, g, -dt);
-
-    /* Find the force. */
-
-    v_scl(F, A, m / dt);
-
-    /* Find the position of the pendulum. */
-
-    v_scl(r, up->E[1], -up->r);
-
-    /* Find the torque on the pendulum. */
-
-    if (fabsf(v_dot(r, F)) > 0.0f)
-        v_crs(T, F, r);
-
-    /* Apply the torque and dampen the angular velocity. */
-
-    v_mad(up->W, up->W, T, dt);
-    v_scl(up->W, up->W,    kd);
-
-    /* Apply the angular velocity to the pendulum basis. */
-
-    sol_rotate(up->E, up->W, dt);
-
-    /* Apply a torque turning the pendulum toward the ball velocity. */
-
-    v_mad(v, up->v, up->E[1], v_dot(up->v, up->E[1]));
-    v_crs(Y, v, up->E[2]);
-    v_scl(Y, up->E[1], 2 * v_dot(Y, up->E[1]));
-
-    sol_rotate(up->E, Y, dt);
-}
-
-/*---------------------------------------------------------------------------*/
-
-/*
- * Compute the states of all switches after DT seconds have passed.
- */
-static void sol_swch_step(struct s_file *fp, float dt)
-{
-    int xi;
-
-    for (xi = 0; xi < fp->xc; xi++)
-    {
-        struct s_swch *xp = fp->xv + xi;
-
-        volatile float t = xp->t;
-
-        if (t < xp->t0)
-        {
-            xp->t = (t += dt);
-
-            if (t >= xp->t0)
-            {
-                int pi = xp->pi;
-                int pj = xp->pi;
-
-                do  /* Tortoise and hare cycle traverser. */
-                {
-                    fp->pv[pi].f = xp->f0;
-                    fp->pv[pj].f = xp->f0;
-
-                    cmd.type        = CMD_PATH_FLAG;
-                    cmd.pathflag.pi = pi;
-                    cmd.pathflag.f  = fp->pv[pi].f;
-                    sol_cmd_enq(&cmd);
-
-                    pi = fp->pv[pi].pi;
-                    pj = fp->pv[pj].pi;
-                    pj = fp->pv[pj].pi;
-                }
-                while (pi != pj);
-
-                xp->f = xp->f0;
-
-                cmd.type          = CMD_SWCH_TOGGLE;
-                cmd.swchtoggle.xi = xi;
-                sol_cmd_enq(&cmd);
-            }
-        }
-    }
-}
-
-/*
- * Compute the positions of all bodies after DT seconds have passed.
- */
-static void sol_body_step(struct s_file *fp, float dt)
-{
-    int i;
-
-    for (i = 0; i < fp->bc; i++)
-    {
-        struct s_body *bp = fp->bv + i;
-        struct s_path *pp = fp->pv + bp->pi;
-
-        volatile float t = bp->t;
-
-        if (bp->pi >= 0 && pp->f)
-        {
-            bp->t = (t += dt);
-
-            if (t >= pp->t)
-            {
-                bp->t  = 0;
-                bp->pi = pp->pi;
-
-                cmd.type        = CMD_BODY_TIME;
-                cmd.bodytime.bi = i;
-                cmd.bodytime.t  = bp->t;
-                sol_cmd_enq(&cmd);
-
-                cmd.type        = CMD_BODY_PATH;
-                cmd.bodypath.bi = i;
-                cmd.bodypath.pi = bp->pi;
-                sol_cmd_enq(&cmd);
-            }
-        }
-    }
-}
-
-/*
- * Compute the positions of all balls after DT seconds have passed.
- */
-static void sol_ball_step(struct s_file *fp, float dt)
-{
-    int i;
-
-    for (i = 0; i < fp->uc; i++)
-    {
-        struct s_ball *up = fp->uv + i;
-
-        v_mad(up->p, up->p, up->v, dt);
-
-        sol_rotate(up->e, up->w, dt);
-    }
-}
-
-/*---------------------------------------------------------------------------*/
-
-static float sol_test_vert(float dt,
-                           float T[3],
-                           const struct s_ball *up,
-                           const struct s_vert *vp,
-                           const float o[3],
-                           const float w[3])
-{
-    return v_vert(T, o, vp->p, w, up->p, up->v, up->r);
-}
-
-static float sol_test_edge(float dt,
-                           float T[3],
-                           const struct s_ball *up,
-                           const struct s_file *fp,
-                           const struct s_edge *ep,
-                           const float o[3],
-                           const float w[3])
-{
-    float q[3];
-    float u[3];
-
-    v_cpy(q, fp->vv[ep->vi].p);
-    v_sub(u, fp->vv[ep->vj].p,
-          fp->vv[ep->vi].p);
-
-    return v_edge(T, o, q, u, w, up->p, up->v, up->r);
-}
-
-static float sol_test_side(float dt,
-                           float T[3],
-                           const struct s_ball *up,
-                           const struct s_file *fp,
-                           const struct s_lump *lp,
-                           const struct s_side *sp,
-                           const float o[3],
-                           const float w[3])
-{
-    float t = v_side(T, o, w, sp->n, sp->d, up->p, up->v, up->r);
-    int i;
-
-    if (t < dt)
-        for (i = 0; i < lp->sc; i++)
-        {
-            const struct s_side *sq = fp->sv + fp->iv[lp->s0 + i];
-
-            if (sp != sq &&
-                v_dot(T, sq->n) -
-                v_dot(o, sq->n) -
-                v_dot(w, sq->n) * t > sq->d)
-                return LARGE;
-        }
-    return t;
-}
-
-/*---------------------------------------------------------------------------*/
-
-static int sol_test_fore(float dt,
-                         const struct s_ball *up,
-                         const struct s_side *sp,
-                         const float o[3],
-                         const float w[3])
-{
-    float q[3];
-
-    /* If the ball is not behind the plane, the test passes. */
-
-    v_sub(q, up->p, o);
-
-    if (v_dot(q, sp->n) - sp->d + up->r >= 0)
-        return 1;
-
-    /* If it's not behind the plane after DT seconds, the test passes. */
-
-    v_mad(q, q, up->v, dt);
-
-    if (v_dot(q, sp->n) - sp->d + up->r >= 0)
-        return 1;
-
-    /* Else, test fails. */
-
-    return 0;
-}
-
-static int sol_test_back(float dt,
-                         const struct s_ball *up,
-                         const struct s_side *sp,
-                         const float o[3],
-                         const float w[3])
-{
-    float q[3];
-
-    /* If the ball is not in front of the plane, the test passes. */
-
-    v_sub(q, up->p, o);
-
-    if (v_dot(q, sp->n) - sp->d - up->r <= 0)
-        return 1;
-
-    /* If it's not in front of the plane after DT seconds, the test passes. */
-
-    v_mad(q, q, up->v, dt);
-
-    if (v_dot(q, sp->n) - sp->d - up->r <= 0)
-        return 1;
-
-    /* Else, test fails. */
-
-    return 0;
-}
-
-/*---------------------------------------------------------------------------*/
-
-static float sol_test_lump(float dt,
-                           float T[3],
-                           const struct s_ball *up,
-                           const struct s_file *fp,
-                           const struct s_lump *lp,
-                           const float o[3],
-                           const float w[3])
-{
-    float U[3] = { 0.0f, 0.0f, 0.0f };
-    float u, t = dt;
-    int i;
-
-    /* Short circuit a non-solid lump. */
-
-    if (lp->fl & L_DETAIL) return t;
-
-    /* Test all verts */
-
-    if (up->r > 0.0f)
-        for (i = 0; i < lp->vc; i++)
-        {
-            const struct s_vert *vp = fp->vv + fp->iv[lp->v0 + i];
-
-            if ((u = sol_test_vert(t, U, up, vp, o, w)) < t)
-            {
-                v_cpy(T, U);
-                t = u;
-            }
-        }
-
-    /* Test all edges */
-
-    if (up->r > 0.0f)
-        for (i = 0; i < lp->ec; i++)
-        {
-            const struct s_edge *ep = fp->ev + fp->iv[lp->e0 + i];
-
-            if ((u = sol_test_edge(t, U, up, fp, ep, o, w)) < t)
-            {
-                v_cpy(T, U);
-                t = u;
-            }
-        }
-
-    /* Test all sides */
-
-    for (i = 0; i < lp->sc; i++)
-    {
-        const struct s_side *sp = fp->sv + fp->iv[lp->s0 + i];
-
-        if ((u = sol_test_side(t, U, up, fp, lp, sp, o, w)) < t)
-        {
-            v_cpy(T, U);
-            t = u;
-        }
-    }
-    return t;
-}
-
-static float sol_test_node(float dt,
-                           float T[3],
-                           const struct s_ball *up,
-                           const struct s_file *fp,
-                           const struct s_node *np,
-                           const float o[3],
-                           const float w[3])
-{
-    float U[3], u, t = dt;
-    int i;
-
-    /* Test all lumps */
-
-    for (i = 0; i < np->lc; i++)
-    {
-        const struct s_lump *lp = fp->lv + np->l0 + i;
-
-        if ((u = sol_test_lump(t, U, up, fp, lp, o, w)) < t)
-        {
-            v_cpy(T, U);
-            t = u;
-        }
-    }
-
-    /* Test in front of this node */
-
-    if (np->ni >= 0 && sol_test_fore(t, up, fp->sv + np->si, o, w))
-    {
-        const struct s_node *nq = fp->nv + np->ni;
-
-        if ((u = sol_test_node(t, U, up, fp, nq, o, w)) < t)
-        {
-            v_cpy(T, U);
-            t = u;
-        }
-    }
-
-    /* Test behind this node */
-
-    if (np->nj >= 0 && sol_test_back(t, up, fp->sv + np->si, o, w))
-    {
-        const struct s_node *nq = fp->nv + np->nj;
-
-        if ((u = sol_test_node(t, U, up, fp, nq, o, w)) < t)
-        {
-            v_cpy(T, U);
-            t = u;
-        }
-    }
-
-    return t;
-}
-
-static float sol_test_body(float dt,
-                           float T[3], float V[3],
-                           const struct s_ball *up,
-                           const struct s_file *fp,
-                           const struct s_body *bp)
-{
-    float U[3], O[3], W[3], u, t = dt;
-
-    const struct s_node *np = fp->nv + bp->ni;
-
-    sol_body_p(O, fp, bp->pi, bp->t);
-    sol_body_v(W, fp, bp->pi, bp->t, dt);
-
-    if ((u = sol_test_node(t, U, up, fp, np, O, W)) < t)
-    {
-        v_cpy(T, U);
-        v_cpy(V, W);
-        t = u;
-    }
-    return t;
-}
-
-static float sol_test_file(float dt,
-                           float T[3], float V[3],
-                           const struct s_ball *up,
-                           const struct s_file *fp)
-{
-    float U[3], W[3], u, t = dt;
-    int i;
-
-    for (i = 0; i < fp->bc; i++)
-    {
-        const struct s_body *bp = fp->bv + i;
-
-        if ((u = sol_test_body(t, U, W, up, fp, bp)) < t)
-        {
-            v_cpy(T, U);
-            v_cpy(V, W);
-            t = u;
-        }
-    }
-    return t;
-}
-
-/*---------------------------------------------------------------------------*/
-
-/*
- * Step the physics forward DT  seconds under the influence of gravity
- * vector G.  If the ball gets pinched between two moving solids, this
- * loop might not terminate.  It  is better to do something physically
- * impossible than  to lock up the game.   So, if we make  more than C
- * iterations, punt it.
- */
-
-float sol_step(struct s_file *fp, const float *g, float dt, int ui, int *m)
-{
-    float P[3], V[3], v[3], r[3], a[3], d, e, nt, b = 0.0f, tt = dt;
-    int c = 16;
-
-    defer_cmds = 1;
-
-    if (ui < fp->uc)
-    {
-        struct s_ball *up = fp->uv + ui;
-
-        /* If the ball is in contact with a surface, apply friction. */
-
-        v_cpy(a, up->v);
-        v_cpy(v, up->v);
-        v_cpy(up->v, g);
-
-        if (m && sol_test_file(tt, P, V, up, fp) < 0.0005f)
-        {
-            v_cpy(up->v, v);
-            v_sub(r, P, up->p);
-
-            if ((d = v_dot(r, g) / (v_len(r) * v_len(g))) > 0.999f)
-            {
-                if ((e = (v_len(up->v) - dt)) > 0.0f)
-                {
-                    /* Scale the linear velocity. */
-
-                    v_nrm(up->v, up->v);
-                    v_scl(up->v, up->v, e);
-
-                    /* Scale the angular velocity. */
-
-                    v_sub(v, V, up->v);
-                    v_crs(up->w, v, r);
-                    v_scl(up->w, up->w, -1.0f / (up->r * up->r));
-                }
-                else
-                {
-                    /* Friction has brought the ball to a stop. */
-
-                    up->v[0] = 0.0f;
-                    up->v[1] = 0.0f;
-                    up->v[2] = 0.0f;
-
-                    (*m)++;
-                }
-            }
-            else v_mad(up->v, v, g, tt);
-        }
-        else v_mad(up->v, v, g, tt);
-
-        /* Test for collision. */
-
-        while (c > 0 && tt > 0 && tt > (nt = sol_test_file(tt, P, V, up, fp)))
-        {
-            cmd.type       = CMD_STEP_SIMULATION;
-            cmd.stepsim.dt = nt;
-            sol_cmd_enq(&cmd);
-
-            sol_body_step(fp, nt);
-            sol_swch_step(fp, nt);
-            sol_ball_step(fp, nt);
-
-            tt -= nt;
-
-            if (b < (d = sol_bounce(up, P, V, nt)))
-                b = d;
-
-            c--;
-        }
-
-        cmd.type       = CMD_STEP_SIMULATION;
-        cmd.stepsim.dt = tt;
-        sol_cmd_enq(&cmd);
-
-        sol_body_step(fp, tt);
-        sol_swch_step(fp, tt);
-        sol_ball_step(fp, tt);
-
-        /* Apply the ball's accelleration to the pendulum. */
-
-        v_sub(a, up->v, a);
-
-        sol_pendulum(up, a, g, dt);
-    }
-
-    sol_cmd_enq_deferred();
-    defer_cmds = 0;
-
-    return b;
-}
-
-/*---------------------------------------------------------------------------*/
-
-int sol_item_test(struct s_file *fp, float *p, float item_r)
-{
-    const float *ball_p = fp->uv->p;
-    const float  ball_r = fp->uv->r;
-
-    int hi;
-
-    for (hi = 0; hi < fp->hc; hi++)
-    {
-        float r[3];
-
-        r[0] = ball_p[0] - fp->hv[hi].p[0];
-        r[1] = ball_p[1] - fp->hv[hi].p[1];
-        r[2] = ball_p[2] - fp->hv[hi].p[2];
-
-        if (fp->hv[hi].t != ITEM_NONE && v_len(r) < ball_r + item_r)
-        {
-            p[0] = fp->hv[hi].p[0];
-            p[1] = fp->hv[hi].p[1];
-            p[2] = fp->hv[hi].p[2];
-
-            return hi;
-        }
-    }
-    return -1;
-}
-
-struct s_goal *sol_goal_test(struct s_file *fp, float *p, int ui)
-{
-    const float *ball_p = fp->uv[ui].p;
-    const float  ball_r = fp->uv[ui].r;
-    int zi;
-
-    for (zi = 0; zi < fp->zc; zi++)
-    {
-        float r[3];
-
-        r[0] = ball_p[0] - fp->zv[zi].p[0];
-        r[1] = ball_p[2] - fp->zv[zi].p[2];
-        r[2] = 0;
-
-        if (v_len(r) < fp->zv[zi].r - ball_r &&
-            ball_p[1] > fp->zv[zi].p[1] &&
-            ball_p[1] < fp->zv[zi].p[1] + GOAL_HEIGHT / 2)
-        {
-            p[0] = fp->zv[zi].p[0];
-            p[1] = fp->zv[zi].p[1];
-            p[2] = fp->zv[zi].p[2];
-
-            return &fp->zv[zi];
-        }
-    }
-    return NULL;
-}
-
-/*
- * Test if the  ball UI is inside a  jump. Return 1 if yes  and fill P
- * with the destination position, return 0 if not, and return 2 if the
- * ball is on the border of a jump.
- */
-int sol_jump_test(struct s_file *fp, float *p, int ui)
-{
-    const float *ball_p = fp->uv[ui].p;
-    const float  ball_r = fp->uv[ui].r;
-    int ji;
-    float l;
-    int res = 0;
-
-    for (ji = 0; ji < fp->jc; ji++)
-    {
-        float r[3];
-
-        r[0] = ball_p[0] - fp->jv[ji].p[0];
-        r[1] = ball_p[2] - fp->jv[ji].p[2];
-        r[2] = 0;
-
-        l = v_len(r) - fp->jv[ji].r;
-        if (l < 0 &&
-            ball_p[1] > fp->jv[ji].p[1] &&
-            ball_p[1] < fp->jv[ji].p[1] + JUMP_HEIGHT / 2)
-        {
-            if (l < - ball_r )
-            {
-                p[0] = fp->jv[ji].q[0] + (ball_p[0] - fp->jv[ji].p[0]);
-                p[1] = fp->jv[ji].q[1] + (ball_p[1] - fp->jv[ji].p[1]);
-                p[2] = fp->jv[ji].q[2] + (ball_p[2] - fp->jv[ji].p[2]);
-
-                return 1;
-            }
-            else
-                res = 2;
-        }
-    }
-    return res;
-}
-
-/*
- * Test and process the event the ball UI enters a switch. Return 1 if
- * a visible  switch is  activated, return 0  otherwise (no  switch is
- * activated or only invisible switches).
- */
-int sol_swch_test(struct s_file *fp, int ui)
-{
-    const float *ball_p = fp->uv[ui].p;
-    const float  ball_r = fp->uv[ui].r;
-    int xi;
-    int res = 0;
-
-    for (xi = 0; xi < fp->xc; xi++)
-    {
-        struct s_swch *xp = fp->xv + xi;
-
-        /* FIXME enter/exit events don't work for timed switches */
-
-        if (xp->t0 == 0 || xp->f == xp->f0)
-        {
-            float l;
-            float r[3];
-
-            r[0] = ball_p[0] - xp->p[0];
-            r[1] = ball_p[2] - xp->p[2];
-            r[2] = 0;
-
-            l = v_len(r) - xp->r;
-
-            if (l < ball_r &&
-                ball_p[1] > xp->p[1] &&
-                ball_p[1] < xp->p[1] + SWCH_HEIGHT / 2)
-            {
-                if (!xp->e && l < - ball_r)
-                {
-                    int pi = xp->pi;
-                    int pj = xp->pi;
-
-                    /* The ball enters. */
-
-                    if (xp->t0 == 0)
-                    {
-                        xp->e = 1;
-
-                        cmd.type         = CMD_SWCH_ENTER;
-                        cmd.swchenter.xi = xi;
-                        sol_cmd_enq(&cmd);
-                    }
-
-                    /* Toggle the state, update the path. */
-
-                    xp->f = xp->f ? 0 : 1;
-
-                    cmd.type          = CMD_SWCH_TOGGLE;
-                    cmd.swchtoggle.xi = xi;
-                    sol_cmd_enq(&cmd);
-
-                    do  /* Tortoise and hare cycle traverser. */
-                    {
-                        fp->pv[pi].f = xp->f;
-                        fp->pv[pj].f = xp->f;
-
-                        cmd.type        = CMD_PATH_FLAG;
-                        cmd.pathflag.pi = pi;
-                        cmd.pathflag.f  = fp->pv[pi].f;
-                        sol_cmd_enq(&cmd);
-
-                        pi = fp->pv[pi].pi;
-                        pj = fp->pv[pj].pi;
-                        pj = fp->pv[pj].pi;
-                    }
-                    while (pi != pj);
-
-                    /* It toggled to non-default state, start the timer. */
-
-                    if (xp->f != xp->f0)
-                        xp->t = 0.0f;
-
-                    /* If visible, set the result. */
-
-                    if (!xp->i)
-                        res = 1;
-                }
-            }
-
-            /* The ball exits. */
-
-            else if (xp->e)
-            {
-                xp->e = 0;
-
-                cmd.type        = CMD_SWCH_EXIT;
-                cmd.swchexit.xi = xi;
-                sol_cmd_enq(&cmd);
-            }
-        }
-    }
-    return res;
-}
-
-/*---------------------------------------------------------------------------*/
diff --git a/share/solid_phys.h b/share/solid_phys.h
deleted file mode 100644 (file)
index cbd15da..0000000
+++ /dev/null
@@ -1,41 +0,0 @@
-/*
- * Copyright (C) 2003 Robert Kooima
- *
- * NEVERBALL is  free software; you can redistribute  it and/or modify
- * it under the  terms of the GNU General  Public License as published
- * by the Free  Software Foundation; either version 2  of the License,
- * or (at your option) any later version.
- *
- * This program is distributed in the hope that it will be useful, but
- * WITHOUT  ANY  WARRANTY;  without   even  the  implied  warranty  of
- * MERCHANTABILITY or  FITNESS FOR A PARTICULAR PURPOSE.   See the GNU
- * General Public License for more details.
- */
-
-#ifndef SOL_PHYS_H
-#define SOL_PHYS_H
-
-#include "solid.h"
-#include "cmd.h"
-
-/*---------------------------------------------------------------------------*/
-
-void sol_body_p(float p[3], const struct s_file *, int pi, float t);
-
-/*---------------------------------------------------------------------------*/
-
-float sol_step(struct s_file *, const float *, float, int, int *);
-
-int   sol_jump_test(struct s_file *, float *, int);
-int   sol_swch_test(struct s_file *, int);
-
-struct s_goal *sol_goal_test(struct s_file *, float *, int);
-int            sol_item_test(struct s_file *, float *, float);
-
-/*---------------------------------------------------------------------------*/
-
-void sol_cmd_enq_func(void (*enq_fn) (const union cmd *));
-
-/*---------------------------------------------------------------------------*/
-
-#endif
diff --git a/share/solid_sim.h b/share/solid_sim.h
new file mode 100644 (file)
index 0000000..f57b2b3
--- /dev/null
@@ -0,0 +1,29 @@
+/*
+ * Copyright (C) 2003 Robert Kooima
+ *
+ * NEVERBALL is  free software; you can redistribute  it and/or modify
+ * it under the  terms of the GNU General  Public License as published
+ * by the Free  Software Foundation; either version 2  of the License,
+ * or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful, but
+ * WITHOUT  ANY  WARRANTY;  without   even  the  implied  warranty  of
+ * MERCHANTABILITY or  FITNESS FOR A PARTICULAR PURPOSE.   See the GNU
+ * General Public License for more details.
+ */
+
+#ifndef SOL_SIM_H
+#define SOL_SIM_H
+
+#include "solid.h"
+
+/*---------------------------------------------------------------------------*/
+
+void sol_init_sim(struct s_file *);
+void sol_quit_sim(void);
+
+float sol_step(struct s_file *, const float *, float, int, int *);
+
+/*---------------------------------------------------------------------------*/
+
+#endif
diff --git a/share/solid_sim_ode.c b/share/solid_sim_ode.c
new file mode 100644 (file)
index 0000000..ea84ae1
--- /dev/null
@@ -0,0 +1,375 @@
+#include <ode/ode.h>
+#include <stdio.h>
+
+#include "solid.h"
+#include "solid_sim.h"
+#include "solid_all.h"
+#include "solid_cmd.h"
+
+#include "array.h"
+#include "vec3.h"
+#include "cmd.h"
+
+/*---------------------------------------------------------------------------*/
+
+static dWorldID      world;
+static dSpaceID      space;
+static dJointGroupID group;
+
+static Array lumps;
+
+static dBodyID *bodies;
+static dBodyID *balls;
+
+/*---------------------------------------------------------------------------*/
+
+/*
+ * Polygons in the same order as planes, one polygon per plane, each
+ * starting with the number of vertices, followed by the vertex
+ * indices in counter-clockwise order.
+ */
+
+struct convex
+{
+    dReal    *planes;
+    unsigned  nplanes;
+    dReal    *points;
+    unsigned  npoints;
+    unsigned *polygons;
+};
+
+static void make_convex(struct convex *cx, struct s_file *fp, struct s_lump *lp)
+{
+    int i;
+
+    /*
+     * A debug build of ODE (0.11.1 as I write this) will print a
+     * bunch of warnings about our planes not containing the origin
+     * and vertices having the wrong winding.  The latter is triggered
+     * by presuming that the former shall be fulfilled, but the former
+     * basically doesn't matter; collision detection handles
+     * origin-facing planes just fine, provided the vertices on the
+     * plane are still counter-clockwise wrt the plane normal.
+     */
+
+    memset(cx, 0, sizeof (*cx));
+
+    if ((cx->planes = malloc(sizeof (dReal) * 4 * lp->sc)))
+    {
+        for (i = 0; i < lp->sc; i++)
+        {
+            struct s_side *sp = fp->sv + fp->iv[lp->s0 + i];
+
+            cx->planes[i * 4 + 0] = sp->n[0];
+            cx->planes[i * 4 + 1] = sp->n[1];
+            cx->planes[i * 4 + 2] = sp->n[2];
+            cx->planes[i * 4 + 3] = sp->d;
+        }
+
+        cx->nplanes = lp->sc;
+    }
+
+    if ((cx->points = malloc(sizeof (dReal) * 3 * lp->vc)))
+    {
+        for (i = 0; i < lp->vc; i++)
+        {
+            struct s_vert *vp = fp->vv + fp->iv[lp->v0 + i];
+
+            cx->points[i * 3 + 0] = vp->p[0];
+            cx->points[i * 3 + 1] = vp->p[1];
+            cx->points[i * 3 + 2] = vp->p[2];
+        }
+
+        cx->npoints = lp->vc;
+    }
+
+    if ((cx->polygons = malloc(sizeof (unsigned) * lp->fc)))
+        for (i = 0; i < lp->fc; i++)
+            cx->polygons[i] = fp->iv[lp->f0 + i];
+}
+
+void free_convex(struct convex *cx)
+{
+    free(cx->planes);
+    free(cx->points);
+    free(cx->polygons);
+}
+
+/*---------------------------------------------------------------------------*/
+
+#define CAT_OBJECT 0x1
+#define CAT_WORLD  0x2
+
+void sol_init_sim(struct s_file *fp)
+{
+    dGeomID geom;
+
+    int i, j;
+
+    dInitODE();
+
+    world  = dWorldCreate();
+    space  = dHashSpaceCreate(0);
+    group  = dJointGroupCreate(0);
+
+    lumps  = array_new(sizeof (struct convex));
+
+    bodies = calloc(fp->bc, sizeof (dBodyID));
+    balls  = calloc(fp->uc, sizeof (dBodyID));
+
+    for (i = 0; i < fp->bc; i++)
+    {
+        struct s_body *bp = fp->bv + i;
+
+        if (bp->pi >= 0)
+        {
+            bodies[i] = dBodyCreate(world);
+            dBodySetKinematic(bodies[i]);
+        }
+        else
+            bodies[i] = 0;
+
+        for (j = 0; j < bp->lc; j++)
+        {
+            struct s_lump *lp = fp->lv + bp->l0 + j;
+            struct convex *cx;
+
+            if (lp->fl & L_DETAIL)
+                continue;
+
+            cx = array_add(lumps);
+
+            make_convex(cx, fp, lp);
+
+            geom = dCreateConvex(space,
+                                 cx->planes, cx->nplanes,
+                                 cx->points, cx->npoints,
+                                 cx->polygons);
+
+            dGeomSetCategoryBits(geom, CAT_WORLD);
+            dGeomSetCollideBits (geom, CAT_OBJECT);
+
+            if (bodies[i]) dGeomSetBody(geom, bodies[i]);
+        }
+
+        if (bodies[i])
+        {
+            float p[3];
+            sol_body_p(p, fp, bp->pi, 0.0f);
+            dBodySetPosition(bodies[i], p[0], p[1], p[2]);
+        }
+    }
+
+    for (i = 0; i < fp->uc; i++)
+    {
+        balls[i] = dBodyCreate(world);
+
+        geom = dCreateSphere(space, fp->uv[i].r);
+
+        dGeomSetCategoryBits(geom, CAT_OBJECT);
+        dGeomSetCollideBits(geom, CAT_WORLD);
+        dGeomSetData(geom, fp->uv + i);
+        dGeomSetBody(geom, balls[i]);
+
+        dBodySetPosition(balls[i],
+                         fp->uv[i].p[0],
+                         fp->uv[i].p[1],
+                         fp->uv[i].p[2]);
+    }
+}
+
+void sol_quit_sim(void)
+{
+    free(balls);
+    free(bodies);
+
+    while (array_len(lumps))
+    {
+        struct convex *cx = array_get(lumps, array_len(lumps) - 1);
+        free_convex(cx);
+        array_del(lumps);
+    }
+
+    array_free(lumps);
+
+    dJointGroupDestroy(group);
+    dSpaceDestroy(space);
+    dWorldDestroy(world);
+
+    dCloseODE();
+}
+
+/*---------------------------------------------------------------------------*/
+
+static const dVector3 v_null = { 0.0, 0.0, 0.0 };
+
+static const dReal *geomvel(dGeomID g)
+{
+    dBodyID b;
+    return (b = dGeomGetBody(g)) ? dBodyGetLinearVel(b) : v_null;
+}
+
+static const dReal *pointvel(dGeomID g, dVector3 p)
+{
+    dBodyID b;
+
+    if ((b = dGeomGetBody(g)))
+    {
+        static dVector3 v;
+        dBodyGetPointVel(b, p[0], p[1], p[2], v);
+        return v;
+    }
+    else
+        return v_null;
+}
+
+/*
+ * Compute the "energy" of the impact, to determine the sound amplitude.
+ */
+static float bump(dContactGeom *contact)
+{
+    float r[3];
+    v_sub(r, geomvel(contact->g1), pointvel(contact->g2, contact->pos));
+    return fabsf(v_dot(contact->normal, r));
+}
+
+static void spin(dContactGeom *contact)
+{
+    struct s_ball *up;
+    const dReal *p, *v, *w;
+    float r[3], d[3];
+
+    /* Why do I know that the first geom is a sphere? */
+
+    up = dGeomGetData(contact->g1);
+
+    p = dGeomGetPosition(contact->g1);
+    v = geomvel(contact->g1);
+    w = pointvel(contact->g2, contact->pos);
+
+    v_sub(r, p, contact->pos);
+    v_sub(d, v, w);
+
+    /* Find the new angular velocity. */
+
+    v_crs(up->w, d, r);
+    v_scl(up->w, up->w, -1.0f / (up->r * up->r));
+}
+
+static void collide(void *data, dGeomID o1, dGeomID o2)
+{
+    dContact contacts[8];
+    dJointID joint;
+
+    int n, i;
+
+    float *b = (float *) data;
+    float  d;
+
+    if ((n = dCollide(o1, o2, 8, &contacts[0].geom, sizeof (dContact))) > 0)
+    {
+        for (i = 0; i < n; i++)
+        {
+            if (*b < (d = bump(&contacts[i].geom)))
+                *b = d;
+
+            spin(&contacts[i].geom);
+
+            contacts[i].surface.mode       = dContactBounce;
+            contacts[i].surface.mu         = 0;
+            contacts[i].surface.bounce     = 0.718;
+            contacts[i].surface.bounce_vel = 0;
+
+            joint = dJointCreateContact(world, group, contacts + i);
+
+            dJointAttach(joint, dGeomGetBody(o1), dGeomGetBody(o2));
+        }
+    }
+}
+
+/*---------------------------------------------------------------------------*/
+
+static void import_state(struct s_file *fp, float dt)
+{
+    int i;
+
+    for (i = 0; i < fp->uc; i++)
+    {
+        dBodySetPosition(balls[i],
+                         fp->uv[i].p[0],
+                         fp->uv[i].p[1],
+                         fp->uv[i].p[2]);
+
+        dBodySetLinearVel(balls[i],
+                          fp->uv[i].v[0],
+                          fp->uv[i].v[1],
+                          fp->uv[i].v[2]);
+
+        dGeomSphereSetRadius(dBodyGetFirstGeom(balls[i]), fp->uv[i].r);
+    }
+
+    for (i = 0; i < fp->bc; i++)
+    {
+        struct s_body *bp = fp->bv + i;
+        float v[3], p[3];
+
+        if (bodies[i])
+        {
+            sol_body_p(p, fp, bp->pi, bp->t);
+            sol_body_v(v, fp, bp->pi, bp->t, dt);
+
+            dBodySetPosition(bodies[i], p[0], p[1], p[2]);
+            dBodySetLinearVel(bodies[i], v[0], v[1], v[2]);
+        }
+    }
+}
+
+static void export_state(struct s_file *fp)
+{
+    int i;
+
+    for (i = 0; i < fp->uc; i++)
+    {
+        struct s_ball *up = fp->uv + i;
+        const dReal *v;
+
+        v = dBodyGetPosition(balls[i]);
+        v_cpy(up->p, v);
+
+        v = dBodyGetLinearVel(balls[i]);
+        v_cpy(up->v, v);
+    }
+}
+
+/*---------------------------------------------------------------------------*/
+
+float sol_step(struct s_file *fp, const float *g, float dt, int ui, int *m)
+{
+    union cmd cmd;
+    float b = 0.0f;
+
+    /* The simulation is advanced only once, commands don't accumulate. */
+
+    sol_cmd_defer = 0;
+
+    import_state(fp, dt);
+
+    dSpaceCollide(space, &b, collide);
+    dWorldSetGravity(world, g[0], g[1], g[2]);
+    dWorldQuickStep(world, dt);
+    dJointGroupEmpty(group);
+
+    sol_body_step(fp, dt);
+    sol_swch_step(fp, dt);
+    sol_ball_step(fp, dt);
+
+    export_state(fp);
+
+    cmd.type       = CMD_STEP_SIMULATION;
+    cmd.stepsim.dt = dt;
+    sol_cmd_enq(&cmd);
+
+    return b;
+}
+
+/*---------------------------------------------------------------------------*/
diff --git a/share/solid_sim_sol.c b/share/solid_sim_sol.c
new file mode 100644 (file)
index 0000000..b74740b
--- /dev/null
@@ -0,0 +1,639 @@
+/*
+ * Copyright (C) 2003 Robert Kooima
+ *
+ * NEVERBALL is  free software; you can redistribute  it and/or modify
+ * it under the  terms of the GNU General  Public License as published
+ * by the Free  Software Foundation; either version 2  of the License,
+ * or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful, but
+ * WITHOUT  ANY  WARRANTY;  without   even  the  implied  warranty  of
+ * MERCHANTABILITY or  FITNESS FOR A PARTICULAR PURPOSE.   See the GNU
+ * General Public License for more details.
+ */
+
+#include <math.h>
+
+#include "vec3.h"
+#include "common.h"
+
+#include "solid.h"
+#include "solid_sim.h"
+#include "solid_all.h"
+#include "solid_cmd.h"
+
+#define LARGE 1.0e+5f
+#define SMALL 1.0e-3f
+
+/*---------------------------------------------------------------------------*/
+/* Solves (p + v * t) . (p + v * t) == r * r for smallest t.                 */
+
+/*
+ * Given vectors A = P, B = V * t, C = A + B, |C| = r, solve for
+ * smallest t.
+ *
+ * Some useful dot product properties:
+ *
+ * 1) A . A = |A| * |A|
+ * 2) A . (B + C) = A . B + A . C
+ * 3) (r * A) . B = r * (A . B)
+ *
+ * Deriving a quadratic equation:
+ *
+ * C . C = r * r                                     (1)
+ * (A + B) . (A + B) = r * r
+ * A . (A + B) + B . (A + B) = r * r                 (2)
+ * A . A + A . B + B . A + B . B = r * r             (2)
+ * A . A + 2 * (A . B) + B . B = r * r
+ * P . P + 2 * (P . V * t) + (V * t . V * t) = r * r
+ * P . P + 2 * (P . V) * t + (V . V) * t * t = r * r (3)
+ * (V . V) * t * t + 2 * (P . V) * t + P . P - r * r = 0
+ *
+ * This equation is solved using the quadratic formula.
+ */
+
+static float v_sol(const float p[3], const float v[3], float r)
+{
+    float a = v_dot(v, v);
+    float b = v_dot(v, p) * 2.0f;
+    float c = v_dot(p, p) - r * r;
+    float d = b * b - 4.0f * a * c;
+
+/* HACK: This seems to cause failures to detect low-velocity collision
+         Yet, the potential division by zero below seems fine.
+    if (fabsf(a) < SMALL) return LARGE;
+*/
+
+    if      (d < 0.0f) return LARGE;
+    else if (d > 0.0f)
+    {
+        float t0 = 0.5f * (-b - fsqrtf(d)) / a;
+        float t1 = 0.5f * (-b + fsqrtf(d)) / a;
+        float t  = (t0 < t1) ? t0 : t1;
+
+        return (t < 0.0f) ? LARGE : t;
+    }
+    else return -b * 0.5f / a;
+}
+
+/*---------------------------------------------------------------------------*/
+
+/*
+ * Compute the  earliest time  and position of  the intersection  of a
+ * sphere and a vertex.
+ *
+ * The sphere has radius R and moves along vector V from point P.  The
+ * vertex moves  along vector  W from point  Q in a  coordinate system
+ * based at O.
+ */
+static float v_vert(float Q[3],
+                    const float o[3],
+                    const float q[3],
+                    const float w[3],
+                    const float p[3],
+                    const float v[3], float r)
+{
+    float O[3], P[3], V[3];
+    float t = LARGE;
+
+    v_add(O, o, q);
+    v_sub(P, p, O);
+    v_sub(V, v, w);
+
+    if (v_dot(P, V) < 0.0f)
+    {
+        t = v_sol(P, V, r);
+
+        if (t < LARGE)
+            v_mad(Q, O, w, t);
+    }
+    return t;
+}
+
+/*
+ * Compute the  earliest time  and position of  the intersection  of a
+ * sphere and an edge.
+ *
+ * The sphere has radius R and moves along vector V from point P.  The
+ * edge moves along vector W from point Q in a coordinate system based
+ * at O.  The edge extends along the length of vector U.
+ */
+static float v_edge(float Q[3],
+                    const float o[3],
+                    const float q[3],
+                    const float u[3],
+                    const float w[3],
+                    const float p[3],
+                    const float v[3], float r)
+{
+    float d[3], e[3];
+    float P[3], V[3];
+    float du, eu, uu, s, t;
+
+    v_sub(d, p, o);
+    v_sub(d, d, q);
+    v_sub(e, v, w);
+
+    /*
+     * Think projections.  Vectors D, extending from the edge vertex Q
+     * to the sphere,  and E, the relative velocity  of sphere wrt the
+     * edge, are  made orthogonal to  the edge vector U.   Division of
+     * the  dot products  is required  to obtain  the  true projection
+     * ratios since U does not have unit length.
+     */
+
+    du = v_dot(d, u);
+    eu = v_dot(e, u);
+    uu = v_dot(u, u);
+
+    v_mad(P, d, u, -du / uu);
+    v_mad(V, e, u, -eu / uu);
+
+    t = v_sol(P, V, r);
+    s = (du + eu * t) / uu; /* Projection of D + E * t on U. */
+
+    if (0.0f <= t && t < LARGE && 0.0f < s && s < 1.0f)
+    {
+        v_mad(d, o, w, t);
+        v_mad(e, q, u, s);
+        v_add(Q, e, d);
+    }
+    else
+        t = LARGE;
+
+    return t;
+}
+
+/*
+ * Compute  the earliest  time and  position of  the intersection  of a
+ * sphere and a plane.
+ *
+ * The sphere has radius R and moves along vector V from point P.  The
+ * plane  moves  along  vector  W.   The  plane has  normal  N  and  is
+ * positioned at distance D from the origin O along that normal.
+ */
+static float v_side(float Q[3],
+                    const float o[3],
+                    const float w[3],
+                    const float n[3], float d,
+                    const float p[3],
+                    const float v[3], float r)
+{
+    float vn = v_dot(v, n);
+    float wn = v_dot(w, n);
+    float t  = LARGE;
+
+    if (vn - wn <= 0.0f)
+    {
+        float on = v_dot(o, n);
+        float pn = v_dot(p, n);
+
+        float u = (r + d + on - pn) / (vn - wn);
+        float a = (    d + on - pn) / (vn - wn);
+
+        if (0.0f <= u)
+        {
+            t = u;
+
+            v_mad(Q, p, v, +t);
+            v_mad(Q, Q, n, -r);
+        }
+        else if (0.0f <= a)
+        {
+            t = 0;
+
+            v_mad(Q, p, v, +t);
+            v_mad(Q, Q, n, -r);
+        }
+    }
+    return t;
+}
+
+/*---------------------------------------------------------------------------*/
+
+/*
+ * Compute the new  linear and angular velocities of  a bouncing ball.
+ * Q  gives the  position  of the  point  of impact  and  W gives  the
+ * velocity of the object being impacted.
+ */
+static float sol_bounce(struct s_ball *up,
+                        const float q[3],
+                        const float w[3], float dt)
+{
+    float n[3], r[3], d[3], vn, wn;
+    float *p = up->p;
+    float *v = up->v;
+
+    /* Find the normal of the impact. */
+
+    v_sub(r, p, q);
+    v_sub(d, v, w);
+    v_nrm(n, r);
+
+    /* Find the new angular velocity. */
+
+    v_crs(up->w, d, r);
+    v_scl(up->w, up->w, -1.0f / (up->r * up->r));
+
+    /* Find the new linear velocity. */
+
+    vn = v_dot(v, n);
+    wn = v_dot(w, n);
+
+    v_mad(v, v, n, 1.7 * (wn - vn));
+
+    v_mad(p, q, n, up->r);
+
+    /* Return the "energy" of the impact, to determine the sound amplitude. */
+
+    return fabsf(v_dot(n, d));
+}
+
+/*---------------------------------------------------------------------------*/
+
+static float sol_test_vert(float dt,
+                           float T[3],
+                           const struct s_ball *up,
+                           const struct s_vert *vp,
+                           const float o[3],
+                           const float w[3])
+{
+    return v_vert(T, o, vp->p, w, up->p, up->v, up->r);
+}
+
+static float sol_test_edge(float dt,
+                           float T[3],
+                           const struct s_ball *up,
+                           const struct s_file *fp,
+                           const struct s_edge *ep,
+                           const float o[3],
+                           const float w[3])
+{
+    float q[3];
+    float u[3];
+
+    v_cpy(q, fp->vv[ep->vi].p);
+    v_sub(u, fp->vv[ep->vj].p,
+          fp->vv[ep->vi].p);
+
+    return v_edge(T, o, q, u, w, up->p, up->v, up->r);
+}
+
+static float sol_test_side(float dt,
+                           float T[3],
+                           const struct s_ball *up,
+                           const struct s_file *fp,
+                           const struct s_lump *lp,
+                           const struct s_side *sp,
+                           const float o[3],
+                           const float w[3])
+{
+    float t = v_side(T, o, w, sp->n, sp->d, up->p, up->v, up->r);
+    int i;
+
+    if (t < dt)
+        for (i = 0; i < lp->sc; i++)
+        {
+            const struct s_side *sq = fp->sv + fp->iv[lp->s0 + i];
+
+            if (sp != sq &&
+                v_dot(T, sq->n) -
+                v_dot(o, sq->n) -
+                v_dot(w, sq->n) * t > sq->d)
+                return LARGE;
+        }
+    return t;
+}
+
+/*---------------------------------------------------------------------------*/
+
+static int sol_test_fore(float dt,
+                         const struct s_ball *up,
+                         const struct s_side *sp,
+                         const float o[3],
+                         const float w[3])
+{
+    float q[3];
+
+    /* If the ball is not behind the plane, the test passes. */
+
+    v_sub(q, up->p, o);
+
+    if (v_dot(q, sp->n) - sp->d + up->r >= 0)
+        return 1;
+
+    /* If it's not behind the plane after DT seconds, the test passes. */
+
+    v_mad(q, q, up->v, dt);
+
+    if (v_dot(q, sp->n) - sp->d + up->r >= 0)
+        return 1;
+
+    /* Else, test fails. */
+
+    return 0;
+}
+
+static int sol_test_back(float dt,
+                         const struct s_ball *up,
+                         const struct s_side *sp,
+                         const float o[3],
+                         const float w[3])
+{
+    float q[3];
+
+    /* If the ball is not in front of the plane, the test passes. */
+
+    v_sub(q, up->p, o);
+
+    if (v_dot(q, sp->n) - sp->d - up->r <= 0)
+        return 1;
+
+    /* If it's not in front of the plane after DT seconds, the test passes. */
+
+    v_mad(q, q, up->v, dt);
+
+    if (v_dot(q, sp->n) - sp->d - up->r <= 0)
+        return 1;
+
+    /* Else, test fails. */
+
+    return 0;
+}
+
+/*---------------------------------------------------------------------------*/
+
+static float sol_test_lump(float dt,
+                           float T[3],
+                           const struct s_ball *up,
+                           const struct s_file *fp,
+                           const struct s_lump *lp,
+                           const float o[3],
+                           const float w[3])
+{
+    float U[3] = { 0.0f, 0.0f, 0.0f };
+    float u, t = dt;
+    int i;
+
+    /* Short circuit a non-solid lump. */
+
+    if (lp->fl & L_DETAIL) return t;
+
+    /* Test all verts */
+
+    if (up->r > 0.0f)
+        for (i = 0; i < lp->vc; i++)
+        {
+            const struct s_vert *vp = fp->vv + fp->iv[lp->v0 + i];
+
+            if ((u = sol_test_vert(t, U, up, vp, o, w)) < t)
+            {
+                v_cpy(T, U);
+                t = u;
+            }
+        }
+
+    /* Test all edges */
+
+    if (up->r > 0.0f)
+        for (i = 0; i < lp->ec; i++)
+        {
+            const struct s_edge *ep = fp->ev + fp->iv[lp->e0 + i];
+
+            if ((u = sol_test_edge(t, U, up, fp, ep, o, w)) < t)
+            {
+                v_cpy(T, U);
+                t = u;
+            }
+        }
+
+    /* Test all sides */
+
+    for (i = 0; i < lp->sc; i++)
+    {
+        const struct s_side *sp = fp->sv + fp->iv[lp->s0 + i];
+
+        if ((u = sol_test_side(t, U, up, fp, lp, sp, o, w)) < t)
+        {
+            v_cpy(T, U);
+            t = u;
+        }
+    }
+    return t;
+}
+
+static float sol_test_node(float dt,
+                           float T[3],
+                           const struct s_ball *up,
+                           const struct s_file *fp,
+                           const struct s_node *np,
+                           const float o[3],
+                           const float w[3])
+{
+    float U[3], u, t = dt;
+    int i;
+
+    /* Test all lumps */
+
+    for (i = 0; i < np->lc; i++)
+    {
+        const struct s_lump *lp = fp->lv + np->l0 + i;
+
+        if ((u = sol_test_lump(t, U, up, fp, lp, o, w)) < t)
+        {
+            v_cpy(T, U);
+            t = u;
+        }
+    }
+
+    /* Test in front of this node */
+
+    if (np->ni >= 0 && sol_test_fore(t, up, fp->sv + np->si, o, w))
+    {
+        const struct s_node *nq = fp->nv + np->ni;
+
+        if ((u = sol_test_node(t, U, up, fp, nq, o, w)) < t)
+        {
+            v_cpy(T, U);
+            t = u;
+        }
+    }
+
+    /* Test behind this node */
+
+    if (np->nj >= 0 && sol_test_back(t, up, fp->sv + np->si, o, w))
+    {
+        const struct s_node *nq = fp->nv + np->nj;
+
+        if ((u = sol_test_node(t, U, up, fp, nq, o, w)) < t)
+        {
+            v_cpy(T, U);
+            t = u;
+        }
+    }
+
+    return t;
+}
+
+static float sol_test_body(float dt,
+                           float T[3], float V[3],
+                           const struct s_ball *up,
+                           const struct s_file *fp,
+                           const struct s_body *bp)
+{
+    float U[3], O[3], W[3], u, t = dt;
+
+    const struct s_node *np = fp->nv + bp->ni;
+
+    sol_body_p(O, fp, bp->pi, bp->t);
+    sol_body_v(W, fp, bp->pi, bp->t, dt);
+
+    if ((u = sol_test_node(t, U, up, fp, np, O, W)) < t)
+    {
+        v_cpy(T, U);
+        v_cpy(V, W);
+        t = u;
+    }
+    return t;
+}
+
+static float sol_test_file(float dt,
+                           float T[3], float V[3],
+                           const struct s_ball *up,
+                           const struct s_file *fp)
+{
+    float U[3], W[3], u, t = dt;
+    int i;
+
+    for (i = 0; i < fp->bc; i++)
+    {
+        const struct s_body *bp = fp->bv + i;
+
+        if ((u = sol_test_body(t, U, W, up, fp, bp)) < t)
+        {
+            v_cpy(T, U);
+            v_cpy(V, W);
+            t = u;
+        }
+    }
+    return t;
+}
+
+/*---------------------------------------------------------------------------*/
+
+/*
+ * Step the physics forward DT  seconds under the influence of gravity
+ * vector G.  If the ball gets pinched between two moving solids, this
+ * loop might not terminate.  It  is better to do something physically
+ * impossible than  to lock up the game.   So, if we make  more than C
+ * iterations, punt it.
+ */
+
+float sol_step(struct s_file *fp, const float *g, float dt, int ui, int *m)
+{
+    float P[3], V[3], v[3], r[3], a[3], d, e, nt, b = 0.0f, tt = dt;
+    int c = 16;
+
+    union cmd cmd;
+
+    sol_cmd_defer = 1;
+
+    if (ui < fp->uc)
+    {
+        struct s_ball *up = fp->uv + ui;
+
+        /* If the ball is in contact with a surface, apply friction. */
+
+        v_cpy(a, up->v);
+        v_cpy(v, up->v);
+        v_cpy(up->v, g);
+
+        if (m && sol_test_file(tt, P, V, up, fp) < 0.0005f)
+        {
+            v_cpy(up->v, v);
+            v_sub(r, P, up->p);
+
+            if ((d = v_dot(r, g) / (v_len(r) * v_len(g))) > 0.999f)
+            {
+                if ((e = (v_len(up->v) - dt)) > 0.0f)
+                {
+                    /* Scale the linear velocity. */
+
+                    v_nrm(up->v, up->v);
+                    v_scl(up->v, up->v, e);
+
+                    /* Scale the angular velocity. */
+
+                    v_sub(v, V, up->v);
+                    v_crs(up->w, v, r);
+                    v_scl(up->w, up->w, -1.0f / (up->r * up->r));
+                }
+                else
+                {
+                    /* Friction has brought the ball to a stop. */
+
+                    up->v[0] = 0.0f;
+                    up->v[1] = 0.0f;
+                    up->v[2] = 0.0f;
+
+                    (*m)++;
+                }
+            }
+            else v_mad(up->v, v, g, tt);
+        }
+        else v_mad(up->v, v, g, tt);
+
+        /* Test for collision. */
+
+        while (c > 0 && tt > 0 && tt > (nt = sol_test_file(tt, P, V, up, fp)))
+        {
+            cmd.type       = CMD_STEP_SIMULATION;
+            cmd.stepsim.dt = nt;
+            sol_cmd_enq(&cmd);
+
+            sol_body_step(fp, nt);
+            sol_swch_step(fp, nt);
+            sol_ball_step(fp, nt);
+
+            tt -= nt;
+
+            if (b < (d = sol_bounce(up, P, V, nt)))
+                b = d;
+
+            c--;
+        }
+
+        cmd.type       = CMD_STEP_SIMULATION;
+        cmd.stepsim.dt = tt;
+        sol_cmd_enq(&cmd);
+
+        sol_body_step(fp, tt);
+        sol_swch_step(fp, tt);
+        sol_ball_step(fp, tt);
+
+        /* Apply the ball's accelleration to the pendulum. */
+
+        v_sub(a, up->v, a);
+
+        sol_pendulum(up, a, g, dt);
+    }
+
+    sol_cmd_enq_deferred();
+    sol_cmd_defer = 0;
+
+    return b;
+}
+
+/*---------------------------------------------------------------------------*/
+
+void sol_init_sim(struct s_file *fp)
+{
+    return;
+}
+
+void sol_quit_sim(void)
+{
+    return;
+}
+
+/*---------------------------------------------------------------------------*/