Algorithm to extract line-like contour in 3d
Here is a generic approach:


Find a point on the intersection of these two curved planes, i.e., on the intersection of $F(\vec{x})=l$ and $G(\vec{x})=s$.  Call that point $\vec{x}_0$.
Trace along the intersection curve, starting at $\vec{x}_0$, and continuing as far as you can.  Then do the same starting from $\vec{x}_0$ in the opposite direction.
Find a new point $\vec{x}_0$ on the intersection of $F(\vec{x})=l$ and $G(\vec{x})=s$ that's not part of what you've traced so far, and go back to step 2.


To instantiate this, we need to know how to implement each of these steps.

To find a point on the intersection of $F(\vec{x})=l$ and $G(\vec{x})=s$, you can use any method for solving a system of nonlinear equations.  There are many.  If you're not familiar with the subject, one simple approach is to minimize the objective function

$$\Phi(\vec{x})= \|F(\vec{x})-l\|^2 + \|G(\vec{x})-s\|^2$$

using, e.g., gradient descent.

To trace along the intersection curve, you can use the gradient to help you.  Suppose that $\vec{x}_0$ is a point on the intersection curve, and expand the Taylor series for $F$ and $G$.  We get

$$F(\vec{x}) = F(\vec{x}_0) + \nabla F(\vec{x}_0) \cdot (\vec{x}-\vec{x}_0) + \dots \approx l + \nabla F(\vec{x}_0) \cdot (\vec{x}-\vec{x}_0)$$

and similarly

$$G(\vec{x}) \approx s + \nabla G(\vec{x}_0) \cdot (\vec{x}-\vec{x}_0).$$

So, for $\vec{x}$ to be on the intersection curve, we require $ \nabla F(\vec{x}_0) \cdot (\vec{x}-\vec{x}_0) = 0$ and $\nabla G(\vec{x}_0) \cdot (\vec{x}-\vec{x}_0)=0$ (assuming $x$ is very near $x_0$).  Since $\nabla F(\vec{x}_0)$, $\nabla G(\vec{x}_0)$, and $\vec{x}_0$ are known and computable, this gives us two linear equations in three unknowns, so we can solve for $\vec{x}$ using a linear equation solver.   In particular, compute the two vectors $\nabla F(\vec{x}_0)$ and $\nabla G(\vec{x}_0)$, then find a unit vector $\vec{d}$ that is orthogonal to both of them.  Now you can trace out the intersection curve from $\vec{x}_0$ by preceding in the direction $\pm \vec{d}$; in other words, the line $\vec{x}_1=\vec{x}_0 + c \vec{d}$ is approximately on the intersection curve, when $c$ is close to 0.  (If you want, you can generate $\vec{x}_1$ in this way, then "fine-tune" it by using it as the starting point for gradient descent with objective function $\Phi$ to compute $\vec{x}_1$ more precisely and to make sure you stay on the intersection curve.)  This gives you a way to step along the intersection curve, one point at a time.

How to implement step 3?  This might depend on the specific functions $F,G$.  One approach is to modify the objective function $\Phi$ to penalize points that are near the part of the curve you've already found, e.g.,

$$\Phi(\vec{x}) = \|F(\vec{x})-l\|^2 + \|G(\vec{x})-s\|^2 + \sum_i {1 \over \|\vec{x}-\vec{x}_i\|^2},$$

where the $x_i$'s are points you've previously found (or some well-spaced-out subset of them).  Another approach might be to do step 1 (gradient descent on $\Phi$, without the penalty term) starting from many different starting points (e.g., spaced on a coarse grid).  There may be a better way in your specific situation; I don't know.

I expect this should be more efficient than enumerating all grid points and evaluating $F,G$ at every grid point, especially as the spacing of your grid points gets small.  It does require you to be able to compute the gradients of $F,G$, but in your situation that seems to be no problem.