If complete integral of differential equation $$ x (p^2 +q^2) = zp $$ ( p is partial derivative of z with respect to x and q is partial derivative of z with respect to y )
Passing through $x=0$ and $z^2 =4y$ then envelope of this family passing through $x=1$ ,$y=1 $ has
As per https://en.wikipedia.org/wiki/Method_of_characteristics#Fully_nonlinear_case, a generalized form of characteristics is obtained by solving the ODE system (Lagrange-Charpit equations) \begin{align} \frac{dx}{2px-z} =\frac{dy}{2xq} =\frac{dz}{x(p^2+q^2)} =\frac{dz}{pz} =\frac{dp}{-q^2} =\frac{dq}{pq} \end{align} The last two ratios give $p^2+q^2=a^2=const.$, thus $a^2x=pz$.
The third last and last ratios give $\dfrac{dz}z=\dfrac{dq}q$ thus $$z=cq\tag{$1z$}$$ which gives $$x=\dfrac{c}{a^2}pq.\tag{$1x$}$$
Then also $\dfrac{dy}{2x}=\dfrac{dq}{p}=\dfrac{cq\,dq}{x}$ which gives $$y=\dfrac{c}{a^2}q^2+d.\tag{$1y$} $$
At this point the characteristic curve is represented as the image of the circle $p^2+q^2=a^2$. Elimination of the partial derivatives $p,q$ implies $$ (ax)^2+(a(y-d))^2=z^2\tag{2} $$ as the relation between $x,y,z$ along a characteristic curve.
The initial curve $x_0=0$, $y_0=\frac14z_0^2$ resp $s\mapsto(0,\frac14s^2,s)$ is part of the solution surface, thus $$ \dot z_0(s)=p_0(s)\dot x_0(s)+q_0(s)\dot y_0(s) \\ \implies 1=p_0·0+q_0·\frac12z_0~~\text{ or }~~q_0z_0=2\tag{3} $$ and from the original equation $$0·(p_0^2+q_0^2)=p_0z_0$$ which requires $p_0=0$ in the general case, thus $q_0^2=a^2$ and $ca^2=cq_0^2=q_0z_0=2$.
From $x=\dfrac{c}{a^2}pq=\dfrac12c^2pq$ we get $$ 4x^2=z^2(z_0^2-z^2)\tag{$4x$} $$ Also $y_0=c+d$ so that $$ y=\frac{z^2}{2}+y_0-\frac{z_0^2}2=\frac{z^2}{2}-\frac{z_0^2}4\tag{$4y$} $$
Thus to find the values on the characteristic curve through $x=1=y$ one needs to solve the system $$ 4=z_0^2z^2-z^4,\\ 4=2z^2-z_0^2. $$ for $z$.
Here is an alternative solution. Our starting point is the complete integral of the PDE, derived by Lutz Lehmann: $$ z^2=a^2[x^2+(y-b)^2]. \tag{1} $$ The intersection of the surface $(1)$ with the plane $x=0$ is the pair of lines $z=\pm a(y-b)$. Now consider the tangent line to the curve $x=0, z^2=4y$ at the point $(0,y_0,\sqrt{4y_0})\,(y_0>0)$: $$ z-\sqrt{4y_0}=z_y(y_0)(y-y_0)=\frac{1}{\sqrt{y_0}}(y-y_0)\implies z=\frac{1}{\sqrt{y_0}}(y+y_0). \tag{2} $$ Comparing $(2)$ with $z=a(y-b)$, we conclude that the condition for the surface $(1)$ to be tangent to the curve $x=0, z^2=4y$ is$^{(*)}$ $$ (a,b)=\left(\frac{1}{\sqrt{y_0}},-y_0\right) \implies b=-\frac{1}{a^2}. \tag{3} $$ Substituting $(3)$ in $(1)$, we obtain the family $$ z^2=a^2\left[x^2+\left(y+\frac{1}{a^2}\right)^2\right]. \tag{4} $$ To find the envelope of this family, we eliminate $a$ between $(4)$ and the equation $$ \frac{\partial}{\partial a}\left\{z^2-a^2\left[x^2+\left(y+\frac{1}{a^2}\right)^2\right]\right\}=0 $$ $$ \implies -2a\left[x^2+\left(y+\frac{1}{a^2}\right)^2\right] +\frac{4}{a}\left(y+\frac{1}{a^2}\right)=0. \tag{5} $$ It follows from $(4)$ and $(5)$ that $$ -\frac{2z^2}{a}+\frac{4}{a}\left(y+\frac{1}{a^2}\right)=0 \implies \frac{1}{a^2}=\frac{z^2}{2}-y, \tag{6} $$ which, substituted in $(4)$, yields $$ z^2\left(\frac{z^2}{2}-y\right)=x^2+\frac{z^4}{4} \implies z^2\left(\frac{z^2}{4}-y\right)=x^2. \tag{7} $$ Solving $(7)$ for $z$ we obtain $$ z=\pm\sqrt{2\left(y+\sqrt{x^2+y^2}\right)}, \tag{8} $$ which yields $z(1,1)=\pm\sqrt{2+2\sqrt{2}}$.
$^{(*)}$ The same relationship between $a$ and $b$ is obtained if we require that $z=-a(y-b)$ is tangent to the negative branch of the curve $z^2=4y$ $($i.e., $z=-\sqrt{4y}).$
Addendum
Eq. $(7)$ also can be derived from Eqs. $(4x)$ and $(4y)$ of Lutz Lehmann's answer: $$ (4x)\implies z_0^2=\frac{4x^2}{z^2}+z^2\quad\text{and}\quad (4y)\implies z_0^2=2z^2-4y, $$ hence $$ \frac{4x^2}{z^2}+z^2=2z^2-4y \implies z^2\left(\frac{z^2}{4}-y\right)=x^2. $$
