A two-dimensional finite difference model for reactive gas assisted laser cutting of metals is presented. The effects due to heat transfer and fluid dynamics are taken into consideration in the present model by simultaneously solving the conservation equations in the molten metal and gas jet regions. The solutions in the two regions are coupled at the erosion front where a free surface model is assumed. The erosion front is allowed to deform so that laser-material coupling can be determined based on the local incident angle. Oxidation reactions are considered only at the erosion front and the heat flux due to oxidation is calculated assuming near-equilibrium oxidation reactions;The solution strategy assumes that the flow above the workpiece can be approximated by an axisymmetric jet impingement problem with the specified nozzle operating conditions. The stagnation conditions on the wall of the workpiece are then used as inlet boundary conditions for the jet flow within the kerf. This approximation avoids the need to solve the near-axisymmetric flow above the workpiece and the flow within the kerf simultaneously;The numerical algorithm employs an implicit formulation. Time accuracy is achieved by Newton-Raphson iteration between time levels. To avoid the decoupling of pressure from the continuity equation in the incompressible Navier-Stokes equations, an artificial compressibility relation is used in place of the continuity equation. The spatial discretization depends on the governing equations of each region. In particular, a central scheme is used for the incompressible Navier-Stokes equations while an upwind flux-difference scheme is used for the complete Navier-Stokes equations. To suppress any oscillatory behavior near flow discontinuities, the flux differences are limited by Roe's Superbee flux limiter;Numerous test cases covering both compressible and incompressible flows are analyzed and presented as a validation step for the computer code. Results for laser cutting of a 6.35 mm thick low carbon steel plate with an oxygen gas jet are presented to demonstrate the capability of the present model to yield the details of the cutting process which are not directly observable in any experiment;In addition to the numerical modeling effort, experiments are performed in the present study to supplement the development of the model. Of particular interest is the characterization of the melt interface profile.