  Based on $(2712.4\pm 14.3)\times10^{6}$ $\psi(3686)$ events, we investigate
four hadronic decay modes of the $P$-wave charmonium spin-singlet state
$h_c(^1P_1) \to h^+ h^- \pi^0/\eta$ ($h=\pi$ or $K$) via the process
$\psi(3686) \to \pi^{0}h_c$ at BESIII. The $h_c \to \pi^+ \pi^- \pi^0$ decay is
observed with a significance of 9.6$\sigma$ after taking into account
systematic uncertainties. Evidences for $h_c \to K^+ K^- \pi^0$ and $h_c \to
K^+ K^- \eta$ are found with significances of $3.5\sigma$ and $3.3\sigma$,
respectively, after considering the systematic uncertainties. The branching
fractions of these decays are measured to be $\mathcal{B}(h_c \to \pi^+ \pi^-
\pi^0)=(1.36\pm0.16\pm0.14)\times10^{-3}$, $\mathcal{B}(h_c \to K^+ K^-
\pi^0)=(3.26\pm0.84\pm0.36)\times10^{-4}$, and $\mathcal{B}(h_c \to K^+ K^-
\eta)=(3.13\pm1.08\pm0.38)\times10^{-4}$, where the first uncertainties are
statistical and the second are systematic. No significant signal of
$h_c\to\pi^+\pi^-\eta$ is found, and the upper limit of its decay branching
fraction is determined to be $\mathcal{B}(h_c\to\pi^+\pi^-\eta) < 4.0 \times
10^{-4}$ at 90% confidence level.
