I am writing a software in which I need to find the shortest distance from an arbitrary point in 3D space to an implicit surface defined by a set of metaballs. I wanted to achieve this by using the Lagrange multiplier, but I am getting stuck on the formula and I hope to find some help. (Please forgive me if I am writing something in a not so correct way, I am mainly a software developer.)
Ok, so here are my prerequisites: A point is defined by the constant vector $[u, v, w]$. My implicit function is the sum over all metaballs defined by: $$\sum_{i=1}^n \frac{r_i ^2}{(x - a_i)^2 + (y - b_i)^2 + (z - c_i)^2} >= 1$$ with $[x, y, z]$ being a point for which I want to evaluate the function, $[a_i, b_i, c_i]$ being the centre of a metaball and $r_i^2$ being the radius of a metaball $m_i$ and $1$ being the threshold for my implicit function. When the value of the sum is $>1$, the point is inside my implicit volume. If the value is exact $1$, then the point in question lies on the implicit surface. Since we are talking about distances, the radius is allways positive, the divisor is allways positive, as well as the whole sum of implicit function values and every partial sum.
Now for the Lagrange multiplier I have the following formulas.
The squared distance (and minimum distance function) from my point to the surface is: $$[d(x, y, z)]^2 = (x - u)^2 + (y - v)^2 + (z - w)^2$$ I am taking the squared distance, since this should make no difference when finding the minimum.
My constraint is the implicit surface equation above: $$g(x, y, z) = (\sum_{i=1}^n \frac{r_i ^2}{(x - a_i)^2 + (y - b_i)^2 + (z - c_i)^2}) - 1 = 0$$
Putting this together with the Lagrange multiplier: $$L(x, y, z, λ) = (x - u)^2 + (y - v)^2 + (z - w)^2 + λ ((\sum_{i=1}^n \frac{r_i ^2}{(x - a_i)^2 + (y - b_i)^2 + (z - c_i)^2}) - 1)$$
Now I am taking the partial derivatives: $$\frac{∂L}{∂x} = 2(x-u)-λ(\sum_{i=1}^n \frac{2r_i ^2(x-a_i)}{((x - a_i)^2 + (y - b_i)^2 + (z - c_i)^2)^2}) =0,$$ $$\frac{∂L}{∂y} = 2(y-v)-λ(\sum_{i=1}^n \frac{2r_i ^2(y-b_i)}{((x - a_i)^2 + (y - b_i)^2 + (z - c_i)^2)^2}) =0,$$ $$\frac{∂L}{∂z} = 2(z-w)-λ(\sum_{i=1}^n \frac{2r_i ^2(z-c_i)}{((x - a_i)^2 + (y - b_i)^2 + (z - c_i)^2)^2}) =0,$$ $$\frac{∂L}{∂λ} = (\sum_{i=1}^n \frac{r_i ^2}{(x - a_i)^2 + (y - b_i)^2 + (z - c_i)^2}) - 1=0$$
Aaaand now I don't know where to go from here. I have a feeling, that I am missing something to make this all a lot less complicated function wise. Ideally I would get to some form of a system of linear equations to put this all into a solver in my code, but I am just not getting to any well working change in my functions to make this work, without, for example, blowing up the functions by multipliying all bases of the sum with the rest. So if anyone has an idea or a hint for me, that would be great!