I am utilitizing rotation vectors (or SORA rotations if you care to call them that) as a means of splitting 3D rotations into three scalar orthogonal variables which are impervious to gimbal lock. (see SO(3)) These variables are exposed to a least-squares optimization algorithm which requires partial derivatives for each function and variable in the system. My optimized functions are all based on object positions, so the derivative I'm computing based on these variables takes the form of a "delta" transformation, of which one component is a rotation vector. I know I need a function D(a, b) which defines a "delta" rotation vector mapping from the rotation space of rotation vector a to that of rotation vector b. One way to derive it is to convert a and b to rotation quaternions qa and qb, and convert qa-1qb back to a rotation vector. This seems roundabout, but perhaps I'm foolish to dismiss it. Given a (hopefully) simplified D, the current rotation r and an axis-vector a I need to compute the limit: limx→0 D(r, r + ax) / x I'm interested in these limits for a = (1, 0, 0), a = (0, 1, 0) and a = (0, 0, 1). The application area is physics simulation and kinematics solving; I'm trying to find the optimal rotation for an object based on various potentially conflicting goals. Thanks in advance for any help.