Choose the origin as a point where a lot of lines intersect, or as a center of a circle, or as a center of mass. Choose the axis so that a lot of lines are parallel to them.
For the specific problem, probably best to choose $P$ as the origin. You might choose $G$, but it has a complex relationship between $A,B,C$, so you'd have to put conditions on those coordinates.
If you wait to substitute for $G=(g_1,g_2)$, you get:
$$\begin{align}GA^2+GB^2+GC^2+3GP^2 = &(g_1-x_1)^2+(g_2-x_2)^2+\\
&(g_1-y_1)^2+(g_2-y_2)^2+\\
&(g_1-z_1)^2+(g_2-z_2)^2+\\
&3g_1^2 + 3g_2^2\\
=&6g_1^2+(x_1^2+y_1^2+z_1^2)-2g_1(x_1+y_1+z_1) +\\
&6g_2^2+(x_2^2+y_2^2+z_2^2)-2g_2(x_2+y_2+z_2)
\end{align}$$
Only here at the end do you use that $x_1+y_1+z_1=3g_1$ and $x_2+y_2+z_2=3g_2$. Then you see that $6g_i^2-2g_i(x_i+y_i+g_i)=0$.
This actually reveals something interesting, by the way - any equality involving the linear combination of squares of lengths and only "linear combinations" of a fixed set of otherwise independent points, is true for all points if any only if it is true on each coordinate - that is, if you assume the points are all on a line. It generalizes to a few other types of expression, as well - mostly, via the law of cosine.
These particular problems work without coordinates, using vector notation. If you define, for $x=(x_1,x_2)$ and $y=(y_1,y_2)$ that $x\cdot y = x_1y_1+x_2y_2$, then the distance between $x$ and $y$, squared, is $(x-y)\cdot(x-y)=x\cdot x -2x\cdot y +y\cdot y$. This "dot product" is nice so you get good properties, like $(a+b)\cdot c = a\cdot c + b\cdot c$, and $a\cdot b = b\cdot a$.
Then if $P=(0,0)$ and $a,b,c$ are any points, and $g=\frac{a+b+c}{3}$, then your first question is:
$$\begin{align}
|AG|^2+|BG|^2+|CG^2|+3|PG|^2=\,& (g-a)\cdot (g-a) + \\
&(g-b)\cdot (g-b) +\\
&(g-c)\cdot (g-c) + \\
&3g\cdot g\\
=\,&6g\cdot g -2g\cdot (a+b+c) + a\cdot a + b\cdot b + c\cdot c\\
=\,&6g\cdot g - 6g\cdot \frac{a+b+c}{3} + a\cdot a + b\cdot b + c\cdot c\\
=\,&6g\cdot g - 6g\cdot g +a\cdot a + b\cdot b + c\cdot c\\
=\,&|PA|^2+|PB|^2+|PC|^2
\end{align}$$
By the way, this equality is actually deeply related to statistics, and can be generalized to any number of points in any number of dimensions:
$$|A_1G|^2+\cdots+|A_nG|^2 + n|PG|^2 =|A_1P|^2+\cdots +|A_nP|^2$$
where $G$ is the center of gravity of the $n$ points $A_1,A_2,\dots,A_n$.
This indicates that the mean of a bunch of values "minimizes" a certain function. If the right side is $f(P)$, then we see it as:
$$f(P)=|GP|^2+f(G)$$
and $P=G$ is the point where $f$ is minimized. The value of $\frac{1}{n}f(G)$ is called the "variance" of $A_1,\dots,A_n$.