A finite-element approach to the numerical solution of the Takagi–Taupin equations expressed in a weak form is presented and applied to simulate the X-ray reflectivity curves, spatial intensity distributions and focusing properties of bent perfect crystals in symmetric reflection geometry. The proposed framework encompasses a new formulation of the Takagi–Taupin equations, which appears to be promising in terms of robustness and stability and supports the Fresnel propagation of the diffracted waves. The presented method is very flexible and has the potential of dealing with dynamical X-ray or neutron diffraction problems related to crystals of arbitrary shape and deformation. The reference implementation based on the commercial COMSOL Multiphysics software package is available to the relevant user community.