Ruby, 126 bytes
->n{[8,32].product(*[(-n..-t=1).map{|i|i*=i;n%i<1&&n/=i;i}*2+[0]]*3).map{|j|d=2-n%2
k,x,y,z=j
2*d*x+y+k*z==n/d&&t+=k-16}
t==1}
found a place to initialize t=1 and expanded the list of squares into a triplet instead of using q to make additional copies.
Ruby , 129 bytes
We check the cases k=8 and k=32 and check if there are twice as many solutions for k=8 than k=32. This is done by adding k-16 to t every time we find a solution. This is either +16 in the case k=32 or -8 in the case k=8. Overall the number is congruent if t is still zerothe same as its initial value at the end of the function.
->n{t=0 #counter for solutions
q=(-n..-1).map{|i|i*=i;n%i<1&&n/=i #make n square free by dividing by -n^2 to -1^2 as necessary
i}*2+[0] #return an array of squares, duplicate for 1^2 to n^2, and add the case 0
[8,32].product(q,q,q).map{|j| #make a cartesian product of all possible values for k,x,y,z and iterate
d=2-n%2 #d=1 for odd n, 2 for even n
k,x,y,z=j #unpack j. k=8,32. x,y,z are the squared values from q.
2*d*x+y+k*z==n/d&&t+=k-16} #test if the current values of k,x,y,z are a valid solution. If so, adjust t by k-16 as explained above.
t==0} #return true if t is zerothe same as its initial value. otherwise false.
Ruby, 129 bytes
We check the cases k=8 and k=32 and check if there are twice as many solutions for k=8 than k=32. This is done by adding k-16 to t every time we find a solution. This is either +16 in the case k=32 or -8 in the case k=8. Overall the number is congruent if t is still zero at the end of the function.
->n{t=0 #counter for solutions
q=(-n..-1).map{|i|i*=i;n%i<1&&n/=i #make n square free by dividing by -n^2 to -1^2 as necessary
i}*2+[0] #return an array of squares, duplicate for 1^2 to n^2, and add the case 0
[8,32].product(q,q,q).map{|j| #make a cartesian product of all possible values for k,x,y,z and iterate
d=2-n%2 #d=1 for odd n, 2 for even n
k,x,y,z=j #unpack j. k=8,32. x,y,z are the squared values from q.
2*d*x+y+k*z==n/d&&t+=k-16} #test if the current values of k,x,y,z are a valid solution. If so, adjust t by k-16 as explained above.
t==0} #return true if t is zero. otherwise false.
Ruby, 126 bytes
->n{[8,32].product(*[(-n..-t=1).map{|i|i*=i;n%i<1&&n/=i;i}*2+[0]]*3).map{|j|d=2-n%2
k,x,y,z=j
2*d*x+y+k*z==n/d&&t+=k-16}
t==1}
found a place to initialize t=1 and expanded the list of squares into a triplet instead of using q to make additional copies.
Ruby , 129 bytes
We check the cases k=8 and k=32 and check if there are twice as many solutions for k=8 than k=32. This is done by adding k-16 to t every time we find a solution. This is either +16 in the case k=32 or -8 in the case k=8. Overall the number is congruent if t is the same as its initial value at the end of the function.
->n{t=0 #counter for solutions
q=(-n..-1).map{|i|i*=i;n%i<1&&n/=i #make n square free by dividing by -n^2 to -1^2 as necessary
i}*2+[0] #return an array of squares, duplicate for 1^2 to n^2, and add the case 0
[8,32].product(q,q,q).map{|j| #make a cartesian product of all possible values for k,x,y,z and iterate
d=2-n%2 #d=1 for odd n, 2 for even n
k,x,y,z=j #unpack j. k=8,32. x,y,z are the squared values from q.
2*d*x+y+k*z==n/d&&t+=k-16} #test if the current values of k,x,y,z are a valid solution. If so, adjust t by k-16 as explained above.
t==0} #return true if t is the same as its initial value. otherwise false.
Ruby, 129 bytes
->n{t=0
[8,32].product(q=(-n..-1).map{|i|i*=i;n%i<1&&n/=i;i}*2+[0],q,q).map{|j|d=2-n%2
k,x,y,z=j
2*d*x+y+k*z==n/d&&t+=k-16}
t==0}
Uses Tunnell's theorem like the other answers. I use a single equation as follows.
2*d*x^2 + y^2 + k*z^2 == n/d where d=2 for even n and d=1 for odd n
We check the cases k=8 and k=32 and check if there are twice as many solutions for k=8 than k=32. This is done by adding k-16 to t every time we find a solution. This is either +16 in the case k=32 or -8 in the case k=8. Overall the number is congruent if t is still zero at the end of the function.
It is necessary to find all solutions to the test equation. I see many answers testing between +/-sqrt n. It is perfectly OK to test also outside these limits if it makes code shorter, but no solutions will be found because the left side of the equation will exceed n. The thing I missed in the beginning is that negative and positive x,y,z are considered separately. Thus -3,0,3 yields three squares 9,0,9 and all solutions must be counted separately (0 must be counted once and 9 must be counted twice.)
Ungolfed code
->n{t=0 #counter for solutions
q=(-n..-1).map{|i|i*=i;n%i<1&&n/=i #make n square free by dividing by -n^2 to -1^2 as necessary
i}*2+[0] #return an array of squares, duplicate for 1^2 to n^2, and add the case 0
[8,32].product(q,q,q).map{|j| #make a cartesian product of all possible values for k,x,y,z and iterate
d=2-n%2 #d=1 for odd n, 2 for even n
k,x,y,z=j #unpack j. k=8,32. x,y,z are the squared values from q.
2*d*x+y+k*z==n/d&&t+=k-16} #test if the current values of k,x,y,z are a valid solution. If so, adjust t by k-16 as explained above.
t==0} #return true if t is zero. otherwise false.