In the context of articulated robotic manipulators, the Forward Kinematics (FK) is a highly non-linear function that maps joint configurations of the robot to poses of its endeffector. Furthermore, while in the most useful cases these functions are neither injective (one-to-one) nor surjective (onto), depending on the robot configuration -- i.e. the sequence of prismatic versus revolute joints, and the number of Degrees of Freedom (DoF) -- the associated Inverse Kinematics (IK) problem may be practically or even theoretically impossible to be solved analytically. Therefore, in the past decades, several approximate methods have been developed for many instances of IK problems. The approximate methods can be divided into two distinct categories: data-driven and numerical approaches. In the first case, data-driven approaches have been successfully used for small workspace domains (e.g., task-driven applications), but not fully explored for large ones, i.e. in task-independent applications where a more general IK is required. Similarly, and despite many successful implementations over the years, numerical solutions may fail if an improper matrix inverse is employed (e.g., Moore-Penrose generalized inverse). In this research, we propose a systematic, robust and accurate numerical solution for the IK problem using the Unit-Consistent (UC) and the Mixed (MX) Inverse methods to invert the Jacobians derived from the Denavit-Hartenberg (D-H) representation of the FK for any robot. As we demonstrate, this approach is robust to whether the system is underdetermined (less than 6 DoF) or overdetermined (more than 6 DoF). We compare the proposed numerical solution to data driven solutions using different robots -- with DoF varying from 3 to 7. We conclude that numerical solutions are easier to implement, faster, and more accurate than most data-driven approaches in the literature, specially for large workspaces as in task-independent applications. We particularly compared the proposed numerical approach against two data-driven approaches: Multi-Layer Perceptron (MLP) and Adaptive Neuro-Fuzzy Inference System (ANFIS), while exploring various architectures of these Neural Networks (NN): i.e. number of inputs, number of outputs, depth, and number of nodes in the hidden layers.