# tracing contour lines of 2d functions

I have a 2-d function $f(x,y)$ for which I can compute the gradient numerically.

The gradient $f_x(x,y){\hat x} + f_y(x,y){\hat y}$ is perpendicular to the contour lines.

From the gradient the unit tangent vector to the contour lines can be computed as $$-\frac{f_y}{(f_x^2 + f_y^2)^{0.5}} {\hat x} + \frac{f_x}{(f_x^2 + f_y^2)^{0.5}} {\hat y}$$

I am trying to compute the contour lines using this unit tangent vector computed above.

The way I proceed is as follows:

1> First get one point such that $f(x_0,y_0) = c$

2> Then take the step as $x_1 = x_0 - t \frac{f_y}{(f_x^2 + f_y^2)^{0.5}}$ and $y_1 = y_0 + t \frac{f_x}{(f_x^2 + f_y^2)^{0.5}}$

(my choice of t is t = 0.0001. The positive value of t traces out the "positive branch" of the contour)

3> I choose t = -0.0001 and repeat the step above which traces out the "negative branch" of the contour.

4> I repeat step 2 and 3 above each time computing the unit tangent vector and taking the steps.

What I realized is keeping t=constant does not help me with the optimal stepping.

My question is:

1> is the method I am using to plot the contour lines correct.

2> how can I dynamically adjust "t" to take the most optimal step and smoothly trace the contour lines.

• You can use adaptive step integrator like RK45. In my experience, it works quite well in tracing out the contour. You may need to correct the trajectory from time to time to ensure it is on the same level set as the initial point, by doing a search in normal direction. Also, the steps taken by RK45 may be too optimal that it skips too much, making your contour not smooth (although the trajectory is correct). Commented Mar 10, 2021 at 14:57

Gradient $$\nabla f(x,y) = \frac{\partial f(x, y)}{\partial x} \hat{e}_x + \frac{\partial f(x, y)}{\partial y} \hat{e}_y = f_x (x,y) \hat{e}_x + f_y (x,y) \hat{e}_y$$ rotated 90° counterclockwise is $$\vec{t}(x, y) = - f_y (x,y) \hat{e}_x + f_x (x,y) \hat{e}_y$$ which is tangent to the contour curve at $(x, y)$; you normalize it to unit length to get $$\hat{t}(x, y) = \frac{\vec{t}(x,y)}{\lVert\vec{t}(x,y)\rVert}$$ So, if you start at some point $( x_0 , y_0 )$ on the desired contour curve, follow $\hat{t}$ exactly, and if $f(x, y)$ is not differentiable everywhere, also follow $-\hat{t}$ exactly starting from $( x_0 , y_0 )$, then you do indeed trace that contour curve.

Unfortunately, as you've already noticed, the problem is the "exactly" part.

Because the contour curve is not composed of line segments, the unit tangent vector $\hat{t}(x, y)$ changes at even infinitesimal distances (for example, consider $f(x, y) = x^2 + y^2$, whose contour curves are circles), there is no easy way to define a step size $\lambda$ such that $$f\left(\vec{p}\right) = f\left(\vec{p} + \lambda\hat{t}(\vec{p})\right)$$ or, using OP's notation, $t$ where $$f(x, y) = f\left( x - t \frac{f_y (x , y)}{\sqrt{(f_x(x,y))^2 + (f_y(x,y))^2}} ,\, y + t \frac{f_x (x , y)}{\sqrt{(f_x(x,y))^2 + (f_y(x,y))^2}} \right)$$ simply because the contour curve curves enough that no finite $\lambda$ ($t$) other than 0 is on the contour curve!

Normally ― for example, when generating contour curves on a topographical map ― a grid of samples is used. While a regular rectangular grid is technically mathematically simpler, a triangular grid produces more "visually pleasing" results. (Humans seem to find defects in a rectangular grid easier than in a hexagonal (triangular) grid.) In this case, the gradient is not needed.

However, knowing the gradient function allows one to use gradient descent/ascent within a grid cell or triangle, locating the minimum and maximum values within the cell or triangle. When minimum and maximum values within each cell or triangle are known, and can be relatively efficiently found for new ones, you can start with a much sparser sampling without the risk of missing entire closed contour curves (since they reside completely within an initial cell or triangle); whenever the minimum and maximum spans the desired contour curve value, you simply subdivide the cell or triangle.

• do you think this formalism can be extended to a 3-D function ? For a 3d function I expect the contour to be a surface in the 3-D space. But I am failing to visualize how f(p) = f(p+lambda * t) can trace out a surface when the point p has 3 coordinates (x,y,z). Commented Jun 13, 2017 at 22:17
• @user1612986: The grid approach is easy, at least if you omit the gradient descent part (and thus sample the field at desired precision using a regular rectangular grid). The linked to is example generic C code that walks a 3D grid in a similar fashion. For simplicity, it produces a point cloud; I omitted the code necessary to turn it into a polygon/triangle mesh (or a Bézier patch set, which is much harder). Commented Jun 13, 2017 at 22:30