After a lot of messing around today I curiously observed that $a^2+b^2$ is only divisible by 3 when both $a$ and $b$ contain factors of 3. I am trying to prove it without using modular arithmetic (because that would be way too easy), but finding it very difficult to do so. Is there an easy way to prove this without using modular arithmetic?
I am also interested in a more general statement. Namely, I want to find the values of $Z$ for which $Z \mid \left( a^2+b^2 \right)$ necessarily implies that $Z \mid \gcd \left( a,b\right)$. We know that $Z$ cannot be 5, because $3^2+4^2=5^2$. More generally, if $Z$ is the largest element of a pythagorean triple then the above implication does not hold.
Hint Compute $a^2+b^2$ modulo $3$ when $a,b$ vary through $0,1,2$. In fact, it suffices you look only at $(1,2),(1,1),(1,0)$, by symmetry.
If I recall correctly: let $p$ be an odd prime. Consider the expression $$f(x,y)=Ax^2+Bxy+Cy^2$$ where $A,B,C$ are integers. Let $\Delta=B^2-4AC$, and suppose $\left(\dfrac{\Delta}{p}\right)=-1$. Then $f(x,y)=0\mod p$ implies $x=y=0\mod p$.
P First, $A=0$ or $C=0$ means $\Delta$ is a q.r. modulo $p$; so $A,C\neq 0$. As $p$ is an odd prime, and $p\not\mid A$, $f(x,y)=0\iff 4Af(x,y)=0$, thus $$(2Ax+By)^2-\Delta y^2=0$$ If $y\neq 0$ then we would get $$\left(\frac{2Ax}y+B\right)^2=\Delta\mod p$$ which is impossible. But then $y=0$ gives $x=0$.
In your case, $\Delta=-4\equiv 2$, so the theorem holds since $\mod 3$ every nonzero is $=1$ upon squaring.