**Abstract** : While the involved physical phenomena are nowadays properly understood, sessile droplet evaporation remains an active research area since several experiments are not comprehensively explained. Indeed, recent experimental works have provided very interesting dynamics [1-5], whose characteristics strongly depend on control parameters, such as substrate temperature, ambient pressure and the various thermo-physical properties of substrate, droplet liquid and surrounding gas. There is no general agreement in the scientific community about the root-causes of the observed behavior. Are they hydro-thermal waves, Rayleigh-Bénard-Marangoni instabilities, or even something else? One of the means to unmistakably answer this question would consist in conducting a stability analysis of the fully coupled and transient fluid flow, heat and mass transfer problem [6]. However, to do so, one must first overcome two major difficulties. The first one concerns the computation of the transient base state, whose complexity arises from the strong non-linearities associated with the three competing driving mechanisms: i) Stefan flow draining liquid towards the interface and more specifically to the triple line when the latter is pinned; ii) buoyancy in the liquid and vapor inducing thermo-solutal convection; iii) surface tension gradient along the liquid-gas interface inducing thermo-soluto-capilary convection. The second level of complexity to perform a relevant stability analysis consists in the intrinsic unsteadiness of the sessile drop evaporation process, in which several time scales co-exist. For particular evaporation configurations, the dynamics of perturbations lies in a much shorter time scale than the total evaporation time, so in this case one can separate perturbations from the base state in a frozen time approach. This enabled us to find that the prevalent perturbations of the base state mainly come from the liquid droplet in which they trigger tridimensional flow structures. Therefore, we have then conducted 3D computations restricted to the droplet and its substrate with a one-sided model, in which heat and mass transfer across the liquid-gas interface is modeled in a semi-analytical way. These unsteady computations are performed during the first stage of the evaporation process in which the triple line is pinned to the substrate, so that the droplet diameter is kept constant, meanwhile the contact angle continuously decreases over time. The obtained results are in good agreement against several experiments performed either under terrestrial conditions or reduced gravity (parabolic flights), cf. figure 1. Figure 1 – 3D unsteady computations: thermal and velocity fields over the droplet surface and in a cut plane. These one-sided model computations helped us to understand one of the potential instability scenarios. At the beginning, a torus roll appears in the wedge close to the pinned triple line and drives eventually an inner roll by viscous stress depending on the droplet aspect ratio (height/radius). As their intensity increases over time they split into thermo-capillary cells owing to Rayleigh-Bénard-Marangoni instability. Then, the number of cells evolves in time owing to the fact that the plane wavelength of these cells is closely related to their depth. As in the first stage of droplet evaporation the triple line is pinned the wetting angle consequently decreases, so does the liquid height, changing correspondingly the thermal Bond number (ratio of Rayleigh to Marangoni numbers). As the number of thermo-convective cells changes over time, hydro-thermal waves take place to redistribute the internal energy in the droplet. In the second stage of evaporation the droplet diameter continuously decreases, so does its perimeter, while the wetting angle stays roughly constant. Therefore the height of peripheral cells remains also constant; consequently their number has to decrease in order to match the droplet perimeter. Here again, hydro-thermal waves enter the game to redistribute the internal energy in the droplet. The dark side of this approach is that computing times are still prohibitive to perform any parametric study this way. Therefore we have conducted a dimensional analysis to account for the main physical components of the problem. As usual, one assumes perturbations