It is hard to tell exactly what level of generalization you are aiming at, as the implicit function theorem is extremely general (holding, for example, for mappings between Banach spaces) and it is unclear to me whether you are seeking a result as strong as that.

But if you are not, here are some thoughts:

1. The specific example you gave can be considered as a special case of the situation where you have a real-valued function on a manifold $f:M\to\mathbb{R}$. Implicit function theorem tells you that when $df(x_0) \neq 0$, then in a neighborhood $U$ of $x_0$ you have $\{x\in U: f(x) = f(x_0)\}$ is a codimension-1 submanifold.

   The case where $df(x_0) = 0$ is dealt with by Morse Theory, when the Hessian (which is well-defined at a critical point) is non-singular. In this case the Morse Lemma gives you the canonical form of $f$ near $x_0$. Note however that in the coordinates provide by the Morse Lemma, the level set $\{x\in U: f(x) = f(x_0)\}$ takes the form of a cone, and so with the exception of $\dim(M) = 2$ you cannot generally find a "implicit function" that is smooth.

2. If you are interested in the case with higher dimensions in the codomain, I think a lot of our current understanding traces to Whitney's paper that invented modern singularity theory. Whitney's paper treated singularities of mappings from the plane and itself, and led to a lot of modern developments. I am not an expert in it, maybe someone else will comment, but at the very least you can just look at papers that cite Whitney to get an idea.

It is hard to tell exactly what level of generalization you are aiming at, as the implicit function theorem is extremely general (holding, for example, for mappings between Banach spaces) and it is unclear to me whether you are seeking a result as strong as that.

But if you are not, here are some thoughts:

1. The specific example you gave can be considered as a special case of the situation where you have a real-valued function on a manifold $f:M\to\mathbb{R}$. Implicit function theorem tells you that when $df(x_0) \neq 0$, then in a neighborhood $U$ of $x_0$ you have $\{x\in U: f(x) = f(x_0)\}$ is a codimension-1 submanifold.

   The case where $df(x_0) = 0$ is dealt with by Morse Theory, when the Hessian (which is well-defined at a critical point) is non-singular. In this case the Morse Lemma gives you the canonical form of $f$ near $x_0$. Note however that in the coordinates provided by the Morse Lemma, the level set $\{x\in U: f(x) = f(x_0)\}$ takes the form of a cone, and so with the exception of $\dim(M) = 2$ you cannot generally find a "implicit function" that is smooth.

2. If you are interested in the case with higher dimensions in the codomain, I think a lot of our current understanding traces to Whitney's paper that invented modern singularity theory. Whitney's paper treated singularities of mappings from the plane to itself, and led to a lot of modern developments. I am not an expert in it, maybe someone else will comment, but at the very least you can just look at papers that cite Whitney to get an idea.