Fingering, veining, channelling and focussing of porous fluids are widely observed phenomena in the Earth's interior, driving a range of geo-processes across all scales. While observations suggest fairly localized flow patterns induced by fractures, the classical Darcian model predicts diffusive behaviour that leads to never-ending spreading and delocalization. We here investigate an alternative physical mechanism without the need to involve fractures. Decompaction weakening leads to the formation and propagation of localized flow-pathways in fluid-saturated porous media. We numerically solve the coupled equations using high-resolution 2-D and 3-D numerical modelling to predict non-linear porous flow in a non-linearly viscously deforming matrix. We show that high-porosity channels may be a dynamic and natural outcome of sufficiently resolved hydromechanical coupling and decompaction weakening. We propose an efficient solution strategy that involves an iterative pseudo-transient numerical method to solve the coupled system of equations in a matrix-free fashion. We discuss benefits and limitations of this approach that performs optimally on hardware accelerators such as graphical processing units and is well-suited for supercomputing. We benchmark the pseudo-transient routines against commonly used direct-iterative solving strategies and show convergence towards identical results. Furthermore, we use the fast solver to systematically study in 2-D the high-porosity channel propagation velocity as a function of bulk and shear viscosity ratios and report discrepancy between 2-D and 3-D configurations. We conclude that the fluid-flow rate in the channels is up to three orders of magnitude higher than expected by pure Darcian flow regimes and show that the high-porosity channels occurrence remains with strain rate dependant shear viscosity. We provide both the 2-D MATLAB-based direct-iterative and pseudo-transient routines for full reproducibility of the presented results and suggest our model configuration as a key benchmark case to validate the implementation of hydromechanical coupling in 2-D and 3-D numerical codes. The routines are available from Bitbucket and the Swiss Geocomputing Centre website, and are also supporting information to this paper.