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.
Thanks in advance.
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.