How to calculate the concentration of all species present in a solution with $\mathrm{0.3~M}$ $\mathrm{NaH_2PO_4}$?
The whole truth is that any time that you add any phosphate ion into an aqueous solution, then you will have all four phosphate species. The $\mathrm{NaH_2PO_4}$ solution contains the ions $\mathrm{Na^+}$, $\mathrm{H^+}$, $\mathrm{OH^−}$, $\mathrm{H_2PO_4^{-}}$, $\mathrm{HPO_4^{2-}}$, $\mathrm{PO_4^{3-}}$ and undissociated acid $\mathrm{H_3PO_4}$.
I form seven independent equations in order to specify the unknown concentrations of the seven species present in an aqueous solution of $\mathrm{NaH_2PO_4}$. these equations are:
$$ K_1 = \frac{[\mathrm{H^+}][\mathrm{H_2PO_4^-}]}{[\mathrm{H_3PO_4}]} = 7.6 \times 10^{−3} \tag{1}$$ $$K_2 = \frac{[\mathrm{H^+}][\mathrm{HPO_4^{2-}}]}{[\mathrm{H_2PO_4^-}]} = 6.2 \times 10^{−8} \tag{2}$$ $$K_{3} = \frac{[\mathrm{H^+}][\mathrm{PO_4^{3-}}]}{[\mathrm{HPO_4^{2-}}]} = 2.1 \times 10^{−13} \tag{3}$$ $$K_\mathrm{w} = {[\mathrm{H^+}][\mathrm{OH^-}]} = 1 \times 10^{−14} \label{eq:4}\tag{4} $$ $$C_\mathrm{P} =[\mathrm{Na^+}]= [\mathrm{NaH_2PO_4}]_0 ={0.3}\tag{5}$$ $$C_\mathrm{P} = [\mathrm{PO_4^{3-}}] + [\mathrm{HPO_4^{2-}}] + [\mathrm{H_2PO_4^-}] + [\mathrm{H_3PO_4}] \tag{6}$$ $${[\mathrm{H_3PO_4}] + [\mathrm{H^+}] = [\mathrm{HPO_4^{2-}}] + 2[\mathrm{PO_4^{3-}}] +\frac{K_\mathrm{w}}{[\mathrm{H^+}]}}\tag{7}$$
Neglecting $[\mathrm{Na^+}]$ and equation $\eqref{eq:4}$, then given $C_\mathrm{P} = {0.3}$ and $K_\mathrm{w}=1\times10^{−14}$, how to solve five equations to find five variables using wolframalpha for solving quintic equation?
For simplicity, write $x_0:=[\texttt{H}_3\texttt{PO}_4]$, $x_1:=[\texttt{H}_2\texttt{PO}_4^-]$, $x_2:=[\texttt{H}\texttt{PO}_4^{2-}]$, and $x_3:=[\texttt{PO}_4^{3-}]$ so that $x_i$ is the concentration of the phosphate species with negative-$i$ charge, for $i=0,1,2,3$. Let $y$ denote $[\texttt{H}^+]$. Therefore, $$x_0+x_1+x_2+x_3=C_P\,,$$ $$x_0+y=x_2+2x_3+\frac{K_w}{y}\,,$$ $$K_1=\frac{yx_1}{x_0}\,,$$ $$K_2=\frac{yx_2}{x_1}\,,$$ and $$K_3=\frac{yx_3}{x_2}\,.$$ This gives $x_1=\dfrac{K_1x_0}{y}$, $x_2=\dfrac{K_2x_1}{y}=\dfrac{K_1K_2x_0}{y^2}$, and $x_3=\dfrac{K_3x_2}{y}=\dfrac{K_1K_2K_3x_0}{y^3}$. Hence, we end up with two equations in $x_0$ and $y$: $$x_0+\frac{K_1x_0}{y}+\frac{K_1K_2x_0}{y^2}+\frac{K_1K_2K_3x_0}{y^3}=C_P\tag{*}$$ and $$x_0+y=\frac{K_1K_2x_0}{y^2}+\frac{2K_1K_2K_3x_0}{y^3}+\frac{K_w}{y}\,.\tag{#}$$ From (*), we easily get $$x_0=\frac{C_P}{1+\frac{K_1}{y}+\frac{K_1K_2}{y^2}+\frac{K_1K_2K_3}{y^3}}\,.$$ Plugging this into (#) yields $$\frac{C_P}{1+\frac{K_1}{y}+\frac{K_1K_2}{y^2}+\frac{K_1K_2K_3}{y^3}}\,\left(1-\frac{K_1K_2}{y^2}-\frac{2K_1K_2K_3}{y^3}\right)+y-\frac{K_w}{y}=0\,.$$ That is, $$y^5+(C_P+K_1)y^4+(K_1K_2-K_w)y^3-(K_1K_w+C_PK_1K_2-K_1K_2K_3)y^2-(K_1K_2K_w+2C_PK_1K_2K_3)y-K_1K_2K_3K_w=0\,.$$ I think you need a numerical solver to solve for $y$ from the last equation. Once you know $y$, you can find the values of $x_0$, $x_1$, $x_2$, and $x_3$.
Mathematica tells me that there is only one positive real solution $y$ (or using Descartes's Rule of Signs, you know that there exists only one positive real solution) to the previous equation, which is $y\approx 2.1\times 10^{-5}$. This means $x_0\approx 8.4\times 10^{-4}$, $x_1\approx 0.30$, $x_2\approx 8.6\times 10^{-4}$, and $x_3\approx 8.5\times 10^{-12}$. In a sense, you should expect $x_1\approx C_P=0.30$. After all, the dissociative reactions into other phosphate species are quite weak.