A procedure to obtain accurate solutions for heat conduction problems in composite media is presented. The governing partial differential equation is cast into an integral form by the application of the boundary integral theory. The resulting integral equation is then solved on a generic element of the problem domain in a typical finite-element fashion. The partitioning and discretization of the problem domain is shown to resolve media heterogeneity in a straightforward and efficient manner. Three methods of weighting internodal thermal conductivity are tested for their influence on the accuracy of the numerical solutions. It is shown that the weightings do exert some influence on the accuracy of solution. The geometric mean taken over two adjacent elements appears to be the most accurate; since the numerical results obtained from this methodology predicted with good agreement the closed form solutions of heat transfer problems exposed to different boundary conditions.