The point of gauge invariance in classical electromagnetics is that you can choose the scalar and vector potential (more appropriate is to say you have one four-vector potential in the sense of special relativity, but unfortunately most textbooks first start with the old-fashioned 1+3-dimensional formalism, so I also use this formalism) with some arbitrariness. Physically meaningful, i.e., measurable quantities are the electric and magnetic components of the electromagnetic field, which are related to the potentials (using Heaviside Lorentz units) via
\vec{E}=-\frac{1}{c} \frac{\partial \vec{A}}{\partial t}-\vec{\nabla} \Phi, \quad \vec{B}=\vec{\nabla} \times \vec{A}.
This ansatz solves the two homogeneous Maxwell equations, and vice versa these Maxwell equations also guarantee the existence of such potentials:
\vec{\nabla} \cdot \vec{B}=0, \quad \frac{1}{c} \frac{\partial \vec{B}}{\partial t}+\vec{\nabla} \times \vec{E}=0.
Any other set of potentials, given by
\Phi'=\Phi+\frac{1}{c} \partial \chi, \quad \vec{A}'=\vec{a}-\vec{\nabla} \chi
obviously gives the same em. field (\vec{E},\vec{B}), no matter which scalar field \chi you take.
Now, usually the task is to solve Maxwell's equations for given charge and current distributions and maybe subject to some boundary conditions. This task can be made much easier by "choosing the gauge" appropriately. This is done by imposing one constraint on the potentials. For dynamical problems most convenient is the Lorenz gauge (sometimes still called Lorentz gauge, which is a bit unjust against Lorenz, who found this gauge way earlier than Lorentz)
\frac{1}{c} \frac{\partial \Phi}{\partial t}+\vec{\nabla} \cdot \vec{A}=0.
Particularly in static or stationary problems the Coulomb gauge is most convenient:
\vec{\nabla} \cdot \vec{A}=0.
There is in principle no way to say the one or the other gauge is more natural than the other.
With the introduction of the potentials you can forget about the homogeneous Maxwell equations and plug the definitions of the potentials into the inhomogeneous Maxwell equations
\vec{\nabla} \cdot \vec{E}=\rho, \quad \vec{\nabla} \times \vec{B} - \frac{1}{c} \frac{\partial \vec{E}}{\partial t}=\vec{j}.<br />
For an arbitrary gauge you get
<br />
-\frac{1}{c} \frac{\partial \vec{\nabla} \cdot \vec{A}}{\partial t} - \Delta \Phi=\rho, \quad<br />
\vec{\nabla} \times (\vec{\nabla} \times \vec{A})+\frac{1}{c^2} \frac{\partial^2 \vec{A}}{\partial t^2}+\vec{\nabla} \frac{1}{c} \frac{\partial \Phi}{\partial t}=\vec{j}.
Of course, the resulting equations change with the choice of gauge, depending on the chosen gauge condition, and this adapts both the vector and the scalar potential to this choice. For the Lorenz gauge the gauge condition leads to
\left (\frac{1}{c^2} \partial_t^2-\Delta \right ) \Phi=\rho, \quad \left (\frac{1}{c^2} \partial_t^2-\Delta \right ) \vec{A}=\vec{j}.
You should figure out, how the equations look in the Coulomb gauge.
It's a good exercise is also to think about the \chi which leads from Lorenz to Coulomb gauge and vice versa. A very nice didactical paper about this subject by J.D. Jackson (the "textbook Jackson") can be found here:
J.D. Jackson, From Lorenz to Coulomb and other explicit gauge transformations, Am. J. Phys. <b>70</b>, 917 (2002)
http://arxiv.org/abs/physics/0204034
