The Acoustic Black Hole (ABH) phenomenon can be capitalized to manipulate and mitigate flexural waves in thin-walled structures. It features unique space-dependent wavenumber variation and wave celerity reduction in the tapered ABH area, thus posing great challenges to the existing modelling methods. In this work, the Partition of Unity Finite Element Method (PUFEM) is revamped to resolve the frequency response of an ABH beam. The method allows incorporating auxiliary interpolation functions in the finite element framework in order to better cope with the ABH oscillating behaviour. Several types of tapered Timoshenko beam elements are constructed by employing enrichment functions based on the ABH wave solutions with the WKB approximations (for general profiles) or the exact solutions (for parabolic profiles). Other enrichment bases, including polynomials, Fourier series and wavelets, are also investigated as hierarchic refinements. Using these enriched elements, structural responses of an ABH beam are computed and compared with the standard FEM. It is shown that the PUFEM can be easily adapted to model ABH effects with a good accuracy and efficiency, outperforming the conventional FEM for solving ABH problems.