Disclaimer: This is answer is given from a mathematical physics point of view, and it is a little bit technical. Any comment or additional answer from other points of view is welcome.
The classical limit of quantum theories and quantum field theories is not straightforward. It is now a very active research topic in mathematical physics and analysis.
The idea is simple: by its own construction, quantum mechanics should reduce to classical mechanics in the limit $\hslash\to 0$. I don't think it is necessary to go into details, however for QM this procedure is now well understood and rigorous form a mathematical standpoint.
For QFTs, such as QED, the situation is similar, although more complicated, and it can be mathematically handled only in few situations. Although it has not been proved yet, I think it is possible to prove convergence to classical dynamics for a (simple) model of QED, describing rigid charges interacting with the quantized EM field.
The Hilbert space is $\mathscr{H}=L^2(\mathbb{R}^{3})\otimes \Gamma_s(\mathbb{C}^2\otimes L^2(\mathbb{R}^3))$ ($\Gamma_s$ is the symmetric Fock space). The Hamiltonian describes an extended charge (with charge/mass ratio $1$) coupled with a quantized EM field in the Coulomb gauge:
\begin{equation*}
\hat{H}=(\hat{p} - c^{-1} \hat{A}(\hat{x}))^2+\sum_{\lambda=1,2}\hslash\int dk\;\omega(k)a^*(k,\lambda)a(k,\lambda)\; ,
\end{equation*}
where $\hat{p}=-i\sqrt{\hslash}\nabla$ and $\hat{x}=i\sqrt{\hslash}x$ are the momentum and position operators of the particle, $a^{\#}(k,
\lambda)$ are the annihilation/creation operators of the EM field (in the two polarizations) and $\hat{A}(x)$ is the quantized vector potential
\begin{equation*}
\hat{A}(x)=\sum_{\lambda=1,2}\int \frac{dk}{(2\pi)^{-3/2}}\;c\sqrt{\hslash/2\lvert k\rvert}\;e_\lambda(k)\chi(k)(a(k,\lambda)e^{ik\cdot x}+a^*(k,\lambda)e^{-ik\cdot x})\; ;
\end{equation*}
with $e_\lambda(k)$ orthonormal vectors such that $k\cdot e_\lambda(k)=0$ (they implement the Coulomb gauge) and $\chi$ is the Fourier transform of the charge distribution of the particle.
The magnetic field operator is $\hat{B}(x)=\nabla\times \hat{A}(x)$ and the (perpendicular) electric field is
$$\hat{E}(x)=\sum_{\lambda=1,2}\int \frac{dk}{(2\pi)^{-3/2}}\;\sqrt{\hslash\lvert k\rvert/2}\;e_\lambda(k)\chi(k)i(a(k,\lambda)e^{ik\cdot x}-a^*(k,\lambda)e^{-ik\cdot x})\; $$
$\hat{H}$ is a self adjoint operator on $\mathscr{H}$, if $\chi(k)/\sqrt{\lvert k\rvert}\in L^2(\mathbb{R}^3)$, so there is a well defined quantum dynamics $U(t)=e^{-it\hat{H}/\hslash}$. Consider now the $\hslash$-dependent coherent states
\begin{equation*}
\lvert C_\hslash(\xi,\pi,\alpha_1,\alpha_2)\rangle=\exp\Bigl(i\hslash^{-1/2}(\pi\cdot x+i\xi \cdot\nabla)\Bigr)\otimes\exp\Bigl(\hslash^{-1/2}\sum_{\lambda=1,2}(a^*_\lambda(\alpha_\lambda)-a_\lambda(\bar{\alpha}_\lambda))\Bigr)\Omega\; ,
\end{equation*}
where $\Omega=\Omega_1\otimes\Omega_2$ with $\Omega_1\in C_0^\infty(\mathbb{R}^3)$ (or in general regular enough, and with norm one) and $\Omega_2$ the Fock space vacuum.
What it should be at least possible to prove is that ($\alpha_\lambda$ is the classical correspondent of $\sqrt{\hslash}a_\lambda$, and it appears inside $E(t,x)$ and $B(t,x)$ below):
\begin{gather*}
\lim_{\hslash\to 0}\langle C_\hslash(\xi,\pi,\alpha_1,\alpha_2),U^*(t)\hat{p}U(t)C_\hslash(\xi,\pi,\alpha_1,\alpha_2)\rangle=\pi(t)\\
\lim_{\hslash\to 0}\langle C_\hslash(\xi,\pi,\alpha_1,\alpha_2),U^*(t)\hat{x}U(t)C_\hslash(\xi,\pi,\alpha_1,\alpha_2)\rangle=\xi(t)\\
\lim_{\hslash\to 0}\langle C_\hslash(\xi,\pi,\alpha_1,\alpha_2),U^*(t)\hat{E}(x)U(t)C_\hslash(\xi,\pi,\alpha_1,\alpha_2)\rangle=E(t,x)\\
\lim_{\hslash\to 0}\langle C_\hslash(\xi,\pi,\alpha_1,\alpha_2),U^*(t)\hat{B}(x)U(t)C_\hslash(\xi,\pi,\alpha_1,\alpha_2)\rangle=B(t,x)\; ;
\end{gather*}
where $(\pi(t),\xi(t),E(t,x),B(t,x))$ is the solution of the classical equation of motion of a rigid charge coupled to the electromagnetic field:
\begin{equation*}
% \left\{
\begin{aligned}
&\left\{\begin{aligned}
\partial_t &B + \nabla\times E=0\\
\partial_t &E - \nabla\times B=-j
\end{aligned}\right. \mspace{20mu}
\left\{\begin{aligned}
\nabla\cdot &E=\rho\\
\nabla\cdot &B=0
\end{aligned}\right.\\
&\left\{\begin{aligned}
\dot{\xi}&= 2\pi\\
\dot{\pi}&= \frac{1}{2}[(\check{\chi}*E)(\xi)+2\pi\times(\check{\chi}*B)(\xi)]
\end{aligned}\right.
\end{aligned}
% \right.
\end{equation*}
with $j=2\pi\check{\chi}(\xi-x)$, and $\rho=\check{\chi}(\xi-x)$ (charge density and current).
To sum up: the time evolved quantum observables averaged over $\hslash$-dependent coherent states converge in the limit $\hslash\to 0$ to the corresponding classical quantities, evolved by the classical dynamics.
Hoping this is not too technical, this picture gives a precise idea of the correspondence between the classical and quantum dynamics for an EM field coupled to a charge with extended distribution (point charges cannot be treated mathematically on a completely rigorous level both classically and quantum mechanically).
Maxwellian electrodynamics fails when quantum mechanical phenomena are involved, in the same way that Newtonian mechanics needs to be replaced in that regime by quantum mechanics. Maxwell's equations don't really "fail", as there is still an equivalent version in QM, it's just the mechanics itself that changes.
Let me elaborate on that one for a bit. In Newtonian mechanics, you had a time-dependent position and momentum, $x(t)$ and $p(t)$ for your particle. In quantum mechanics, the dynamical state is transferred to the quantum state $\psi$, whose closest classical analogue is a probability density in phase space in Liouvillian mechanics. There are two different "pictures" in quantum mechanics, which are exactly equivalent.
In the Schrödinger picture, the dynamical evolution is encoded in the quantum state $\psi$, which evolves in time according to the Schrödinger equation. The position and momentum are replaced by static operators $\hat x$ and $\hat p$ which act on $\psi$; this action can be used to find the expected value, and other statistics, of any measurement of position or momentum.
In the Heisenberg picture, the quantum state is fixed, and it is the operators of all the dynamical variables, including position and momentum, that evolve in time, via the Heisenberg equation.
In the simplest version of quantum electrodynamics, and in particular when no relativistic phenomena are involved, Maxwell's equations continue to hold: they are exactly the Heisenberg equations on the electric and magnetic fields, which are now operators on the system's state $\psi$. Thus, you're formally still "using" the Maxwell equations, but this is rather misleading as the mechanics around it is completely different. (Also, you tend to work on the Schrödinger picture, but that's beside the point.)
This regime is used to describe experiments that require the field itself to be quantized, such as Hong-Ou-Mandel interferometry or experiments where the field is measurably entangled with matter. There is also a large gray area of experiments which are usefully described with this formalism but do not actually require a quantized EM field, such as the examples mentioned by Anna. (Thus, for example, black-body radiation can be explained equally well with discrete energy levels on the emitters rather than the radiation.)
This regime was, until recently, pretty much confined to optical physics, so it wasn't really something an electrical engineer would need to worry about. That has begun to change with the introduction of circuit QED, which is the study of superconducting circuits which exhibit quantum behaviour. This is an exciting new research field and it's one of our best bets for building a quantum computer (or, depending on who you ask, the model used by the one quantum computer that's already built. ish.), so it's something to look at if you're looking at career options ;).
The really crazy stuff comes in when you push electrodynamics into regimes which are both quantum and relativistic (where "relativistic" means that the frequency $\nu$ of the EM radiation is bigger than $c^2/h$ times the mass of all relevant material particles). Here quantum mechanics also changes, and becomes what's known as quantum field theory, and this introduces a number of different phenomena. In particular, the number of particles may change over time, so you can put a photon in a box and come back to find an electron and a positron (which wouldn't happen in classical EM).
Again, here the problem is not EM itself, but rather the mechanics around it. QFT is built around a concept called the action, which completely determines the dynamics. You can also build classical mechanics around the action, and the action for quantum electrodynamics is formally identical to that of classical electrodynamics.
This regime includes pair creation and annihilation phenomena, and also things like photon-photon scattering, which do seem at odds with classical EM. For example, you can produce two gamma-ray beams and make them intersect, and they will scatter off each other slightly. This is inconsistent with the superposition principle of classical EM, as it breaks linearity, so you could say that the Maxwell equations have failed - but, as I pointed out, it's a bit more subtle than that.
Best Answer
Here's a stab at a short and intuitive answer (i.e. math lite). Warning, I ignore basically all factors that aren't equal to 1 so sorry if I get something wrong. $\hbar = \pi = 1/2 = 2 = 4 = 1...$ etc.
Solutions to Maxwell's equations can be decomposed into monochromatic electromagnetic waves*. An electromagnetic wave has an (complex) amplitude $a\propto E$ which oscillates in time. It can be shown that the energy (or energy flux depending on the exact physical situation) in the electromagnetic wave is proportional to the squared amplitude of the wave:
$$ H = \omega a^*a $$
Instead of expressing the amplitude as a single complex number we can express it in terms of the real an imaginary parts of that number:
$$ a = \text{Re}(a) + i \text{Im}(a) = x + i p $$
$$ H/\omega = x^2 + p^2 $$
WERE this a mass on a spring harmonic oscillator we would identify $x$ and $p$ as the canonically conjugate position and momentum variables for the classical harmonic oscillator. These are known to have a Poission bracket like
$$ \{x, p\} = 1 $$
Canonical quantization works by elevating the canonically conjugate variables to quantum operators (can be thought of as operators on a Hilbert space) with commutation relations:
$$ [\hat{x}, \hat{p}] = i $$
It then follows that
$$ [\hat{a}, \hat{a}^{\dagger}] = 1 $$
from which one can derive that
$$ H = \omega\hat{a}^{\dagger} \hat{a} = \omega\hat{n} $$
Where the eigenvalues of $\hat{n}$ are non-negative integers. This is what is meant by the energy of the harmonic oscillator being quantized. Note that this implies that the complex amplitude of the harmonic oscillator is quantized as well.
Back to the electromagnetic field. It turns out** that 1) electromagnetism can be cast in a Lagrangian form, thus giving us a context in which we can discuss canonically conjugate variables and 2) The real and imaginary parts of the electromagnetic field amplitude, $a$, given by $x$ and $p$, can in fact be related to a pair of conjugate variables of electromagnetism. These conjugate variables are related to the electric field and the time derivative of the electric field. There are also relationships, through Maxwell's equations, to the magnetic field.
This means that we are justified in making the theory of electromagnetism quantum by letting $a\rightarrow \hat{a}$ with
$$ [\hat{a},\hat{a}^{\dagger}]=1 $$
This leads us to the conclusion that the energy in a single mode of the electromagnetic is quantized in units of $\omega$, and likewise the amplitude of the electromagnetic wave is also quantized.
This is what a photon is. A photon is a quantized excitation in the electomagnetic wave, in the same way we can have a quantized excitation in the quantum harmonic oscillator. Also, at this point I'll point out that it is possible to have excitations in the classical electromagnetic field as well. If you ever find yourself feeling funny about photons try to spend more time thinking about the similarities, rather than the differences, between classical excitations of the electromagnetic field (or any oscillator) and quantum excitations of the same.
As described in the comments, the photon arises from electromagnetism through an application of canonical quantization to the canonically conjugate variables which make up a single-mode monochromatic solution to Maxwell's equations.
*These are often taken to be plane waves, but, depending on the boundary conditions for the EM problem, we can decompose the general solution to Maxwell's equations into many different families of spatial modes. Any such decomposition suffices for the description of a photon. The shape of the photon is given by the spatial pattern of the spatial mode that we are quantizing. That is photons can have different shapes, and a photon in one mode can be decomposed as a superposition of photons in other modes. There, I answered about 100 physics stack exchange questions about the shape of a photon in this footnote!
**My favorite references for this are Quantum and Atom Optics by Daniel Steck and UC Berkeley Physics 221A/B Lecture Notes by Robert Littlejohn.