The double spherical harmonics angular approximation in the lowest order, i.e. double P{sub 0} (DP{sub 0}), is developed for the solution of time-dependent non-equilibrium grey radiative transfer problems in planar geometry. The standard P{sub 1} angular approximation represents the angular dependence of the radiation specific intensity using a linear function in the angular domain -1 {le} {mu} {le} 1. In contrast, the DP{sub 0} angular approximation represents the angular dependence as isotropic in each half angular range -1 {le} {mu} < 0 and 0 < {mu} {le} 1. Neglecting the time derivative of the radiation flux, both the P{sub 1} and DP{sub 0} equations can be written as a single diffusion equation for the radiation energy density. Although the DP{sub 0} diffusion approximation is expected to be less accurate than the P{sub 1} diffusion approximation at and near thermodynamic equilibrium, the DP{sub 0} angular approximation can more accurately capture the complicated angular dependence near the non-equilibrium wave front. We develop an adaptive angular technique that locally uses either the DP{sub 0} or the P{sub 1} diffusion approximation depending on the degree to which the radiation and material fields are in thermodynamic equilibrium. Numerical results are presented for a test problem due to Su and Olson for which a semi-analytic transport solution exists. The numerical results demonstrate that the adaptive P{sub 1}-DP{sub 0} diffusion approximation can yield improvements in accuracy over the standard P{sub 1} diffusion approximation for non-equilibrium grey radiative transfer
