Localized Boundary-Domain Singular Integral Equations Based on Harmonic Parametrix for Divergence-Form Elliptic PDEs with Variable Matrix Coefficients