The work deals with a stable finite element discretisation of Dirichlet problem for extended Brinkman-Forchheimer equations describing flows in fixed bed reactors in two and three dimensional domains. Since the considered model problem is rare, it has not been intensively studied by other authors. Its nonlinear character and nonstandard form of differential operators occurring in the equations require a special treatment. The existence and uniqueness of weak solutions are established in suitable Sobolev spaces. Stability and a priori estimates for finite element pairs with discontinuous pressure approximation are presented for low Reynolds numbers. In the case of high Reynolds numbers equal order approximation and local projection stabilisation method are used. The enhancement of accuracy of finite element solutions on axis parallel grids is carried out using superconvergence phenomena. In the last part results concerning discrete maximum principle are established for the low order local projection scheme applied to scalar convection-diffusion-reaction problems. The presented numerical results are in a good agreement with the developed theory.