  Inclusive $\chi_{cJ}$ $(J=0,1,2)$ production from $\Upsilon(1S)$ decay is
studied within the framework of nonrelativistic QCD (NRQCD) factorization at
leading order in $v_Q^2$, which includes the contributions of
$b\bar{b}({}^3S_1^{[1]})\to c\bar{c}(^3P_J^{[1]})+X$ and
$b\bar{b}({}^3S_1^{[1]})\to c\bar{c}(^3S_1^{[8]})+X$. For both channels, the
short-distance coefficients are calculated through ${\cal O}(\alpha_s^5)$,
which is next-to-leading order for the second one. By fitting to the measured
$\Upsilon(1S)$ branching fractions to $\chi_{c1}$ and $\chi_{c2}$, we obtain
the color-octet long-distance matrix element (LDME)
$\langle\mathcal{O}^{\chi_{c0}}({}^3S_1^{[8]})\rangle
=(4.04\pm0.47_{-0.34}^{+0.67})\times10^{-3}$ GeV$^3$, where the first error is
experimental and the second one due to the renormalization scale dependence, if
we use as input $\langle\mathcal{O}^{\chi_{c0}}({}^3P_0^{[1]})\rangle=0.107$
GeV$^5$ as obtained via potential-model analysis. Previous LDME sets, extracted
from data of prompt $\chi_{cJ}$ hadroproduction, yield theoretical predictions
that systematically undershoot or mildly overshoot the experimental values of
$\mathcal{B}(\Upsilon\to \chi_{cJ}+X)$.
