An inverse finite element method is developed for simultaneous solution of multi-dimensional solid-liquid phase boundaries and associated three-dimensional solid phase temperature fields. The technique, applicable to quasisteady phase change problems, fixes element nodes at known temperature locations and uses a coarse, spatially limited mesh. This approach is designed to : (1) reduce direct and overall solution costs, (2) eliminate iterative direct solutions associated with temperature dependent thermophysical properties, (3) limit calculations to the heat affected zone and (4) eliminate ad hoc assumptions concerning the boundary heat flux distribution. The inverse algorithm couples a nonlinear solid phase conduction solver with conjugate gradient minimization. First-order regularization and upwind differencing are implemented to improve solution smoothness and stability and an analog welding experiment is used to investigate the technique’s capabilities.