# How to solve this analytic geometry problem completely

I want to find a plane that passes through points {1,0,0} and {0,1,0} and is tangent to surface $$z(x,y)=x^{2}+y^{2}$$.

Solve[{a, b, c}.{1, 0, 0} == d && a*0 + b*1 + c*0 == d &&
a*x0 + b*y0 + c*z0 == d && z0 == x0^2 + y0^2 &&
VectorAngle[{a, b, c}, {-2 x0, -2 y0, 1}] ==
0,(*MatrixRank[{2x0,2y0,1},{a,b,c}]\[Equal]1*){a, b, c, d, x0, y0,
z0}]


But I can't get the answer I want with the above code(the answer is $$z=0$$ and $$2x+2y-z=2$$). What should I do?

Clear["*"];
f = x^2 + y^2 - z;
Solve[{Grad[f, {x, y, z}].({x, y, z} - {1, 0, 0}) == 0 ,
Grad[f, {x, y, z}].({x, y, z} - {0, 1, 0}) == 0, f == 0}, {x, y, z}]
Grad[f, {x, y, z}].({X, Y, Z} - {x, y, z}) == 0 /. % // Simplify


This is pretty easy to do with the extant region functionality:

eq = {u, v, u^2 + v^2};
plane = InfinitePlane[eq, Transpose[D[eq, {{u, v}}]]];

sols = plane /. Solve[RegionMember[plane, #] & /@ {{1, 0, 0}, {0, 1, 0}}, {u, v}]
{InfinitePlane[{0, 0, 0}, {{1, 0, 0}, {0, 1, 0}}],
InfinitePlane[{1, 1, 2}, {{1, 0, 2}, {0, 1, 2}}]}


If you want to see the plane equations themsleves, you need an extra step:

Simplify[RegionMember[#, {x, y, z}], {x, y, z} ∈ Reals] & /@ sols
{z == 0, 2 + z == 2 (x + y)}


I expect two solutions!

First we consider the plane p0,p1,p2

p1 = {1, 0, 0};
p2 = {0, 1, 0};
p0 = {x0, y0, x0^2 + y0^2};


with normal

en= Cross[p0 - p1, p0 - p2];


The normal is forced to be pependicular to the tangentplane at p0

p0/. Solve[{n.D[p0, x0] == 0, n.D[p0, y0] == 0}, {x0, y0}, Reals]
(*{{0, 0, 0}, {1, 1, 2}}*)

GraphicsRow[{Show[{Plot3D[x^2 + y^2, {x, -3, 3}, {y, -3, 3},Mesh -> None,BoxRatios -> {1, 1, 1}]
, Graphics3D[{Point[{p0, p1, p2}], InfiniteLine[{p1, p2}], InfinitePlane[{p1, p2, p0}]} /. sol[[1]]]}],
Show[{Plot3D[x^2 + y^2, {x, -3, 3}, {y, -3, 3}, Mesh -> None,BoxRatios -> {1, 1, 1}]
, Graphics3D[{Point[{p0, p1, p2}], InfiniteLine[{p1, p2}], InfinitePlane[{p1, p2, p0}]} /. sol[[2]]]}]}]


Try this:

Gxyz = z - x^2 - y^2;
p = {x, y, z};
p1 = {1, 0, 0};
p2 = {0, 1, 0};
equ1 = n.(p - p1) == 0
equ2 = n.(p - p2) == 0
equ3 = Gxyz == 0
sol = Solve[{equ1, equ2, equ3}, p]
n0 = n /. sol

gr1 = ContourPlot3D[Gxyz == 0, {x, -1, 1}, {y, -1, 1}, {z, -1, 1}];
gr2 = ContourPlot3D[n0[[1]].(p - p1) == 0, {x, -1, 1}, {y, -1, 1}, {z, -1, 1}];
gr3 = ContourPlot3D[n0[[2]].(p - p1) == 0, {x, -1, 1}, {y, -1, 1}, {z, -1, 1}];
Show[gr1, gr2, gr3]

Solve[{Resolve[
ForAll[{x, y}, x^2 + y^2 + a*x + b*y + c >= 0] &&
Exists[{x, y}, x^2 + y^2 + a*x + b*y + c == 0],
Reals], {a, b, -1}.{1, 0, 0} + c == 0, {a, b, -1}.{0, 1, 0} + c ==
0}, {a, b, c}]
a*x + b*y - z + c == 0 /. %
`