## 2 Simulation Method
### Steady-State Vacuum Simulation  
The technique to simulate the steady-state evaporation/condensation by maintaining the temperature gradient within the liquid film placed in the periodic box was proposed by [33]. In contrast, we propose a simulation approach to simulate the steady-state evaporation process into virtual vacuum by maintaining a constant number of atoms within the liquid film during the simulations. As in case of deletion vacuum simulation, the periodic boundary conditions are applied only in \(x\) and \(y\) directions. The atoms at the left- and right-side interfaces are managed in the following way (see Figure 4). The right-side interface atoms that cross the vapor boundary are considered as evaporated, and their information (position, speed etc.) is registered. The evaporated atoms are immediately transported to the left-side vapor boundary position (without changing any attribute of the atom except for the position in \(z\) axis), thus, creating the virtual vacuum condition near the right-side liquid surface. Then, the evaporated and transported atoms form the mass flux \(J_{tr}=J_{sp}\), which replenish the liquid phase from the left-side of the liquid film, thus, maintaining the number of atoms within the liquid film constant. The transportation of the evaporated atoms slightly shifts the center of mass of the system after each evaporation; therefore, the system is recentered to the initial center of mass position every simulation timestep without changing any relative positions or velocities between the atoms to compensate for the center of mass shift. The transportation of evaporated atoms creates a problem with possible atom overlapping near the left-side vapor boundary. This problem is solved in the following way. The mirror reflection boundary is placed near the left-side vapor boundary as shown in Figure 4. Thus, the spontaneously evaporated atoms from the left-side liquid surface are reflected to the liquid phase after the collisions with the incoming transported atoms or the reflection boundary, and the volume near the left-side vapor boundary is left empty of the interface atoms. However, the overlapping of the short-range repulsion force occurs every 0.22 ns in the volume beyond the reflection boundary when newly evaporated atoms are transported on top of the previously evaporated atoms. Therefore, the dynamics of the atoms in the limit NVE region is solved in the NVE ensemble with the velocity limit, which is allowed for atoms to reach in this region. The velocity limit dissipates the potential energy to a certain value when it is converted to kinetic energy during the interactions  
Figure 3: Deletion vacuum simulation scheme with denoted interface region boundary positions and mass fluxes through the interfaces. The density profile for the scheme was obtained from the production run of the liquid-vapor phase equilibria simulation. between the overlapped atoms. The velocity limit of 800 m/s set in our simulations affects only the highly overlapped atoms due to transportation, while the dynamics of regular atom collisions in the region is solved in the unaltered NVE ensemble. This way, each overlap in the limit NVE region does not break the simulation but rather results in slight growth of the total kinetic energy, which is later accounted for with the thermostat applied to the liquid region. The normalized and absolute positions for steady-state vacuum simulation boundaries (listed in Figure 4) are given in Table 3. In addition, the spontaneous evaporation rates \(J_{sp}\) in all vacuum simulations were obtained in the same way as the mass flux \(I_{out}\) in the phase equilibria simulation (see Equation (6)).