Constructing higher-order difference schemes are always challenging for boundary value problems. The core part is to define boundary
enclosure in such a way that guarantees stability and uniform order of accuracy for all nodes. In this work, we develop sixth-order implicit finite
difference scheme for 2-D heat conduction equation with Dirichlet boundary conditions. The computed generalized eigenvalues of implicit finite
difference matrices have negative real parts that guarantees stability in the
case of Crank-Nicolson method. The validity of our developed numerical
scheme is clearly reflected by the numerical testing.