How to compute a Differential Equation With Dot Product. I have the following differential equation, that I would like to plug into a numerical integrator.
$$
\frac{d \mathbf r}{dt} \cdot \hat{n}(\mathbf r) = 0
$$
Where $\mathbf r(t)$, is the trajectory of a particle and $\hat n$ is a normal vector to a surface. So this equation says that the velocity of the particle must be perpendicular to the normal of the surface. It would be useful to have $d\mathbf r / dt$ alone. To achieve this we impose other restrictions like:
$$
\left|\left| \frac{d \mathbf r}{dt} \right|\right| = 1
$$
We get two equations:
$$
\begin{pmatrix}
n_x & n_y & n_z\\
u_x & u_y & u_z\\
1 & 1 & 1
\end{pmatrix}
\begin{pmatrix}
u_x \\
u_y \\
u_z \\
\end{pmatrix} =
\begin{pmatrix}
0 \\
1 \\
u_x + u_y + u_z \\
\end{pmatrix}
$$
And we let the last equation be whatever we want, obtaining the following result:
$$
\mathbf u(t) \equiv \frac{d \mathbf r}{dt}
$$
$$
\mathbf u(t) =\begin{pmatrix}
u_z - u_y & n_y - n_z & u_y n_z - u_z n_y\\
u_x - u_z & n_z - n_x & u_z n_x - u_x n_z\\
u_y - u_x & n_x- n_y & u_x n_y - u_y n_x
\end{pmatrix} \begin{pmatrix}
0 \\
1 \\
u_x + u_y + u_z \\
\end{pmatrix} = 
\begin{pmatrix}
n_y -n_z + (u_x + u_y +u_z) (u_yn_z-u_zn_y) \\
n_z -n_x + (u_x + u_y +u_z) (u_yn_z-u_zn_y) \\
n_x -n_y + (u_x + u_y +u_z) (u_yn_z-u_zn_y) \\
\end{pmatrix}
$$
$$
\mathbf u(t) = (\mathbf i + \mathbf j + \mathbf k) \times \hat n + ((\mathbf i + \mathbf j + \mathbf k) \cdot \mathbf u(t)) (\mathbf u(t) \times\hat n) 
\\
(\mathbf i + \mathbf j + \mathbf k) \equiv \mathbf e \\
\mathbf e \cdot \mathbf u(t) \equiv k(t) \\
\mathbf u(t) = (\mathbf e  + k(t)\mathbf u(t)) \times \hat n
\\
$$
This last equation is equivalent (I suppose), to the first one I proposed.
But how would I plug this into a program. What I'm seeking for is an equation of the form,
$$
\mathbf u_{i + 1} = F(\mathbf u_i, \mathbf r_i)
$$
that relates the previous $\mathbf u(t)$ with the new one, but I do not have $\mathbf u(t)$ alone (basically because I can't). How I'm supposed to do this?
 A: Let the surface be defined by $f(x,y,z)=0$, and the path on the surface be $r(t)=(x(t),y(t),z(t))$. I will sketch how to derive the geodesic equations for the components of $r$ in terms of $f$, assuming there are no external forces present. For the full details see here (from where this answer is taken, essentially verbatim). Form the vector
$$\tag{1}
N=\frac{\nabla f}{|\nabla f|}\\
$$
The geodesic curvature is
$$\tag{2}
\kappa_g=\ddot{r}\cdot (N\times \dot{r})
$$
When $r(t)$ is a geodesic, $\kappa_g=0$. Substitute (1) into (2) to find
$$\tag{3}
\varepsilon_{ijk}\ddot{r}_i\dot{r}_k\frac{\partial f}{\partial r_j}=0
$$
Where $\varepsilon$ is the Levi-Civita symbol and summation convention is used. We need two more coupled differential equations. They are
$$\tag{4}
\ddot{r}\cdot \dot{r}=0
$$
And,
$$\tag{5}
\frac{d}{dt}\frac{df}{dt}=0=\text{calculus...}
$$
Equations (3) (4) (5) may be solved for $\ddot{r}_i$ (it is done explicitly in the link above). To put these into a numerical solver, you can write them as six coupled first order differential equations. So let $p_i=\dot{r}_i$ and our equation are in the variables $\{ x,y,z,p_x,p_y,p_z\}$
$$
\left[\matrix{\dot{p_x} \\\dot{p_y}\\ \dot{p_y}\\ \dot{x}\\ \dot{y} \\ \dot{z}}\right]=\left[ \matrix{\ddot{x} \\ \ddot{y} \\ \ddot{z} \\ p_x \\ p_y \\p_z} \right]
$$
Note that $\ddot{x}$, $\ddot{y}$ and $\ddot{z}$ on the RHS are to be replaced by the lengthy expressions which were solved for by using eqs (3) (4) and (5). You will need to specify an initial position $r(0)$ and velocity $\dot{r}(0)$.
A: There is not enough data to compute the trajectory. The equation you present just says that the velocity will always be tangent to the surface where the motion occurs. This must be supplemented with the actual equations of motion. A natural choice would be
$$
m\, r''(t) = (0, 0, -g).
$$
