Desalination has been recognized as a climate-independent source of producing potable water that has a great importance due to the lack of drinking water in many countries. One of the methods extensively used for removing salt from seawater is thermal distillation, which is known as a reliable technology for producing high-quality water through evaporation. Among different evaporation methods, a liquid falling film has been widely developed in the thermal desalination industry. The flow characteristics of the falling evaporation process are complicated, and heat and mass transfer are strongly interdependent. Some previous numerical studies simplified the solution by using a boundary layer approximation, which is not capable of modelling reverse flows and recirculation zones. The previous models that solved the full governing equations (elliptic approach) applied the volume of fluid method (VOF), which is not able to precisely determine the location of the liquid and gas interface, therefore, it cannot accurately apply the interface jump conditions. Additionally, these studies commonly used mass transfer models based on estimations and empirical correlations, which introduce limitations to numerical simulations. Another limitation found in the previous studies is that none of them has considered the occurrence of the supersaturation phenomenon in the mixture phase, which may occur during the evaporation process. In this research, a two-dimensional, two-phase model was developed using an in-house Fortran code that solves the steady-state, full elliptic governing equations for the liquid and gas phases. This modeling represents a significant advancement in the simulation of the evaporation phenomenon that employs innovative methods to address the aforementioned problems. This model is believed to be the most comprehensive and detailed approach compared to any previous numerical studies on falling film evaporation. A variety of boundary conditions and geometries have been utilized to analyze the details of evaporation and desalination phenomena. The results were then compared to available experimental and numerical data for validation purposes. In each stage, a new set of detailed results, including the effects of different parameters variation on flow characteristics and heat and mass transfer during the evaporation process are presented.