mathematica is not able to solve this system

Question

I have this system of three equation that mathematica does not resolve.

x = 1005790/1000000;
y = 3012034/1000000;
varx = 9979050/1000000;
vary = 12985206/1000000;
covxy = 9986968/1000000;
z = 3614929/1000000;
varz = 25466937/1000000;
z3 = 323087800/1000000;
eq1 = z == av + at*x + a1*y 
eq2 = varz == at^2 *varx + a1^2 *vary + 2*covxy *a1*at 
eq3 = (av + x*at + a1*y)^3 +3*(av + x*at + a1*y)*(at^2 *varx + a1^2 * vary + 2* at*a1* covxy) == z3 
FindRoot[{eq1 , eq2, eq3}, {{av, 0.4}, {at, 0.7}, {a1, 1.1}}]
Reduce[{eq1 , eq2, eq3}, {av, at, a1}]

Reduce return False and FindRoot is printing the error: FindRoot::jsing: Encountered a singular Jacobian at the point {av,at,a1} = {0.4,0.7,1.1}. Try perturbing the initial point(s). But even if I change the initial points by a Little, The program is unable to solve the equations. I have plotted the three equations minus the known term, and graphically it seems that there is a solution (I fixed av to the correct value, and I see that the three functions are near zero at the same point).

code to plot:

prim[av_, at_, a1_] := av + at*x + a1*y - z; 
sec [av_, at_, a1_] := at^2 *varx + a1^2 *vary + 2*covxy *a1*at - varz;
terz [av_, at_, a1_] := (av + x*at + a1*y)^3 + 
   3*(av + x*at + a1*y)*(at^2 *varx + a1^2 * vary + 
      2* at*a1* covxy) - z3;
Plot3D[{prim[0.1, at, a1], sec[0.1, at, a1], terz[0.1, at, a1]}, {at, 
  0.4, 0.6}, {a1, 0.9, 1.1}]

Answer 1

Try NMinimize:

eq = {-z + av + at*x + a1*y,-varz + at^2*varx + a1^2*vary + 2*covxy*a1*at, (av + x*at + a1*y)^3 +3*(av+ x*at + a1*y)*(at^2*varx + a1^2*vary +2*at*a1*covxy) - z3};

sol=NMinimize[eq . eq, {a1, at, av}]
(*{8.30607*10^-6, {a1 -> 0.966855, at ->0.539657, av -> 0.157078}}*)
eq /. sol[[2]] 
(* {-0.00286933, -0.000269104, 0.0000248339}*)

addendum

If we look examplary for positive solutions, two ContourPlot's show

sol1 = Solve[eq[[2]] == 0, a1][[-1]]
Column[{
  ContourPlot[eq[[2]] == 0, {at, 0, 2}, {a1, 0, 2}, 
   FrameLabel -> {at, a1}], 
  ContourPlot[
   Evaluate[Thread[eq[[{1, 3}]] == 0 /. sol1]], {at, 0, 2}, {av, 0, 
    2}, ContourStyle -> {Blue, {Dashed, Red}}, 
   FrameLabel -> {at, av}]}]

that solutions a1[at],av [at]exist for ~.5<at<1.5