The following code
Block[
{ifun1, values},
ifun1 = Interpolation[
{{{0.0, 0.0}, 0.0}, {{0.0, 1.0}, 0.1},
{{1.0, 0.0}, 0.2}, {{1.0, 1.0}, 1.0}},
InterpolationOrder -> 1, Method -> "Spline"
];
Print[ContourPlot[ifun1[x, y], {x, 0, 1}, {y, 0, 1}]];
];
produces this contour plot
why are the contours smooth? I would expect Mathematica to pick particular grid triangles in the interpolation and was expecting to see those.

ContourPlot[x y, {x, 0, 1}, {y, 0, 1}]$\endgroup$InterpolationOrder -> 1means the interpolants will be linear in each variable — bilinear in this case — not a linear polynomial. $\endgroup$ifun1 = {{{0.0, 0.0}, 0.0}, {{0.0, 1.0}, 0.1}, {{1.0, 0.0}, 0.2}, {{1.0, 1.0}, 1.0}} // NDSolve`FEM`ElementMeshInterpolation[{NDSolve`FEM`ToElementMesh[#[[All, 1]]]}, #[[All, 2]]] &;? $\endgroup$Method -> "Spline"only works for regular grid, and your test sample happens to form a regular grid. If it's not, for example:Block[{ifun1, values}, ifun1 = Interpolation[{{{10^-6 + 0.0, 0.0}, 0.0}, {{0.0, 1}, 0.1}, {{1.0, 0.0}, 0.2}, {{1.0, 1.0}, 1.0}}]; ContourPlot[ifun1[x, y], {x, 0, 1}, {y, 0, 1}]]Interpolationautomatically turns to the unstructured method and you get an unsmooth result. $\endgroup$