For your last question, Yes!, You can use gradient penalty with LSGAN too and you can read this. Gradient penalty is a trick and independent of what divergence/ distance you use.
First one is longer though. TO reiterate, WGAN is trained with Wasserstein distance, not divergence! This is important as divergence is a weaker notion of distance due to the fact that divergence is not symmetric. Ex. $KL(p || q) \neq KL(q || p)$and in fact exploits different properties!
I'm gonna explain it as in the chronological order to understand why WGAN is important.
Before WGAN, GANs(Vanilla GAN, DCGAN , LSGAN and many other GANs before WGAN) , were trained to minimize an f-divergence (KL, JSD, Pearson...). If we take derivative of JSD with respect to generator parameters while real and generator data distribuitons are far from each other, gradient converges to zero. Very bad generator! Pearson divergence provides gradient to generators even if distributions are far away.
Wasserstein OR Kantorovich-Rubinstein metric OR Earth Mover's Distance is the distance between two continuous probability distributions defined as
$$ W(p_r, p_g) = \inf_{\gamma \sim \Pi(p_r, p_g)} \mathbb{E}_{(x, y) \sim \gamma}[\| x-y \|] $$
where $\Pi(p_r, p_g)$ is the set of all possible joint probability distributions between real and generator data distribution. $\gamma \sim \Pi(p_r, p_g)$ defines
So, what makes Wasserstein different than others?
Referring to WGAN paper, say we have two distributions, $\textit{P}$ and $\textit{Q}$ :
$$ \forall (x, y) \in P, x = 0 \text{ and } y \sim U(0, 1)\\
\forall (x, y) \in Q, x = \theta, 0 \leq \theta \leq 1 \text{ and } y \sim U(0, 1)\\$$
When $\theta \neq 0$ , there's no overlap:
]1)
$$
\begin{aligned}
D_{KL}(P \| Q) &= \sum_{x=0, y \sim U(0, 1)} 1 \cdot \log\frac{1}{0} = +\infty \\
D_{KL}(Q \| P) &= \sum_{x=\theta, y \sim U(0, 1)} 1 \cdot \log\frac{1}{0} = +\infty \\
D_{JS}(P, Q) &= \frac{1}{2}(\sum_{x=0, y \sim U(0, 1)} 1 \cdot \log\frac{1}{1/2} + \sum_{x=0, y \sim U(0, 1)} 1 \cdot \log\frac{1}{1/2}) = \log 2\\
W(P, Q) &= |\theta|
\end{aligned}
$$
Wasserstein provides a smooth measure even if distributions are far from each other. This helps for a stable learning procedure, eliminating mode collapse, and improving the learnable class of manifolds (check this).
However, people don't use the Wasserstein metric as it is due to intractability of computing infimum. Using Kantorovich-Rubinstein duality:
$$
W(p_r, p_g) = \frac{1}{K} \sup_{\| f \|_L \leq K} \mathbb{E}_{x \sim p_r}[f(x)] - \mathbb{E}_{x \sim p_g}[f(x)]
$$
to measure least upper bound on the function. Function needs to be K-Lipschitz continuous(Strongly advising this to read.) .
Skipping some formalities ,say our function is from a family of K-Lipschitz continuois functions , $\{ f_w \}_{w \in W}$ , parameterized by $w$, Wasserstein distance is measured by :
$$
L(p_r, p_g) = W(p_r, p_g) = \max_{w \in W} \mathbb{E}_{x \sim p_r}[f_w(x)] - \mathbb{E}_{z \sim p_r(z)}[f_w(g_\theta(z))]
$$