A middleground is to use fixed timesteps with the closed form solution. Then the integration factors can be statically precomputed as constants, and every frame becomes just a 2x2 matrix multiplication.

And even if you are to calculate the integration factors based on dynamic frame rates, they only need to be computed just once per frame for all objects that share the same damping configuration.