  We present a variational augmentation procedure to optimize the exponents of
Gaussian continuum basis sets for simulating strong-field laser ionization
phenomena such as higher harmonic generation (HHG) in atoms and ions using the
time-dependent configuration interaction (TDCI) method. We report the
distribution of the optimized exponents and discuss how efficiently the
resulting basis functions span the variational space to describe the
near-continuum states involved in HHG. Further, we calculated the higher
harmonic spectra of three two-electron systems -- H$^{-}$, He and Li$^{+}$ --
generated by 800nm driving laser-pulses with pulse-width of 54fs and peak
intensities in the tunnel ionization regime of each system. We analyze the
performance of these basis sets with an increasing number of higher angular
momentum functions and show that up to $g$-type functions are required to
obtain qualitatively accurate harmonic spectra. Additionally, we also comment
on the impact of electron correlation on the HHG spectra. Finally, we show that
by systematically augmenting additional shells we model the strong-field
dynamics at higher laser peak intensities.
