Boundary value methods are applied to find transient solutions of MM2 queueing systems with two heterogeneous servers under a variant vacation policy. An iterative method is employed to solve the resulting large linear system and a Crank-Nicolson preconditioner is used to accelerate the convergence. Numerical results are presented to demonstrate the efficiency of the proposed method.