The prices of assets differ when considering transaction costs and their illiquidity, thus affecting their returns. As volatility derivatives such as variance swaps, gamma swaps, corridor variance swaps and volatility swaps depend on the expected values of these returns, we develop a high-order algorithm so that differences in the fair values of the volatility derivatives under different transaction costs and illiquidity models can be precisely observed. A high-order finite-difference algorithm is used to solve a framework of partial differential equations. A local mesh refinement strategy is applied at the singularities that appear in the pricing problem of some volatility derivatives so that nonuniform five-point stencil finite differences can be implemented in space. The time-stepping is dealt with by an efficient iterative scheme coupled with a Richardson extrapolation and Carathéodory–Fejér approximations. Consequently, the highly convergent solutions are used to analyze the effect of transaction costs and illiquidity on the fair values of volatility derivatives, and the results demonstrate an increase in these values when different transaction costs and liquidity models are implemented.