A method for solving the inverse geometrical problem is presented by reconstructing the unknown subsurface cavity geometry using boundary element methods, a genetic algorithm, and Nelder-Mead non-linear simplex optimization. The heat conduction problem is solved utilizing the boundary element method, which calculates the difference between the measured temperature at the exposed surface and the computed temperature under the current update of the unknown subsurface flaws and cavities. In a first step, clusters of singularities are utilized to solve the inverse problem and to identify the location of the centroid(s) of the subsurface cavity(ies)/flaw(s). In a second step, the reconstruction of the estimated cavity(ies)/flaw(s) geometry(ies) is accomplished by utilizing an anchored grid pattern upon which cubic spline knots are restricted to move in the search for unknown geometry. Solution of the inverse problem is achieved using a genetic algorithm accelerated with the Nelder-Mead non-linear simplex. To optimize the cubic spline interpolated geometry, the flux (Neumann) boundary conditions are minimized using a least squares functional. The automated algorithm successfully reconstructs single and multiple subsurface cavities within two dimensional mediums. The solver is also shown to accurately predict cavity geometries with random noise in the boundary condition measurements. Subsurface cavities can be difficult to detect based on their location. By applying different boundary conditions to the same geometry, more information is supplied at the boundary, and the subsurface cavity is easily detected despite its low heat signature effect at the boundaries. Extensions to three-dimensional applications are outlined