  We propose an analytical approach to high-harmonic generation (HHG) for
nonperturbative low-frequency and high-intensity fields based on the
(Jeffreys-)Wentzel-Kramers-Brillouin (WKB) approximation. By properly taking
into account Stokes phenomena of WKB solutions, we obtain wavefunctions that
systematically include the repetitive dynamics of production and acceleration
of electron-hole pairs and quantum interference due to phase accumulation
between different pair production times (St\"{u}ckelberg phase). Using the
obtained wavefunctions without relying on any phenomenological assumptions, we
explicitly compute electric current (including intra- and inter-band
contributions) as the source of HHG for a massive Dirac system in
(1+1)-dimensions under an ac electric field. We demonstrate that the WKB
approximation agrees well with numerical results obtained by solving the
time-dependent Schr\"{o}dinger equation and point out that the quantum
interference is important in HHG. We also predict in the deep nonperturbative
regime that (1) harmonic intensities oscillate with respect to electric-field
amplitude $E_0$ and frequency $\Omega$, with a period determined by the
St\"{u}ckelberg phase; (2) the cutoff order of HHG is determined by
$2eE_0/\hbar \Omega^2$, with $e$ being the electron charge; and that (3)
non-integer harmonics, controlled by the St\"{u}ckelberg phase, appear as a
transient effect. Our WKB theory is particularly suited for a parameter regime,
where the Keldysh parameter $\gamma=(\Delta/2)\Omega/eE_0$, with $\Delta$ being
the gap size, is small. This parameter regime corresponds to intense lasers in
the terahertz regime for realistic massive Dirac materials. Our analysis
implies that the so-called HHG plateau can be observed at the terahertz
frequency within the current technology.
