Analysis of the Atomization Process by Molecular Dynamics Simulation

Investigation of the atomizer flows has received a lot of attentions due to their broad applications in many fields, e.g. industrial processes, agriculture, meteorology, medicine, etc. Historically, less attention has been devoted to the internal flow development in atomizers because of the small size of practical atomizers, which makes measurements difficult. However, internal flows in atomizers are of interest due to their potential effect on the subsequent atomization process that takes place external to the atomizers.


Introduction
Investigation of the atomizer flows has received a lot of attentions due to their broad applications in many fields, e.g. industrial processes, agriculture, meteorology, medicine, etc. Historically, less attention has been devoted to the internal flow development in atomizers because of the small size of practical atomizers, which makes measurements difficult. However, internal flows in atomizers are of interest due to their potential effect on the subsequent atomization process that takes place external to the atomizers.
Owing to the difficulties inherent in the simulation of the atomizer internal and external flows, no systematic numerical researches have been carried out on this subject. One way to study internal flows in atomizers is to conduct experimental investigations by the use of large-scale models. Another approach is to do detailed numerical simulations of the flow. An advantage of the numerical simulation is that practical-scale atomizers can be studied as easily as large scale ones. The methods of analyzing atomizer flows can be divided into two categories, namely, the macroscopic analysis and the microscopic analysis, each with its own advantages and disadvantages. The former approach is usually adopted as an aided tool for atomizer design due to its relatively lower computational cost. However, when atomization occurs, the liquid breaks into ultra-fine nano-scale droplets. Under such situation, the Navier-Stokes equation based on continuum concept is not suitable and the microscopic analysis should be used for subsequent simulation of the atomization process. In the past decade, the author has carried out a series of researches on the macroscopic analysis of the atomizer internal and external flows (Yeh, 2002(Yeh, , 2003(Yeh, , 2004(Yeh, , 2005(Yeh, , 2007. The macroscopic behavior and characteristics of the atomizer flows, such as discharge coefficient, spray angle, film thickness, atomizer geometry, etc. were analyzed. However, as mentioned above, the macroscopic analysis is not suitable when atomization occurs. The author's previous studies were therefore initially confined to the investigation of the early stage of the atomization process, during which the liquid had not yet ruptured to ultra-fine nano-scale droplets. In order to improve this inherent deficiency of the macroscopic analysis, the author has employed molecular dynamics (MD) simulation to investigate the microscopic evolution of the atomizer internal and external flows. In the past few years, the author has studied the atomization processes of a nano-scale liquid thread and a nanojet. Liquid thread and nanojet are two of the most fundamental and important phenomena during atomization process. Although the size of a nano-scale atomizer is much smaller than 12 6 (r) 4 rr where r denotes the distance between two molecules, and σ are the representative scales of energy and length, respectively. The Lennard-Jones fluid in this research is taken to be argon for its ease of physical understanding. The parameters for argon are as follows (Kawano, 1998) : the length parameter σ=0.354 nm, the energy parameter /k B =93.3K, and the molecular weight m=6.64×10 -26 kg, where k B =1.38×10 -23 J/K denotes the Boltzmann constant. The cut-off radius r c beyond which the intermolecular interaction is neglected is 5.0σ.
The simulation domain is schematically shown in Fig.1, with periodic boundary conditions applied in all three directions. Simulation domain dimensions, temperatures and number of molecules, are listed in Table 1, together with some simulation results. The time integration of motion is performed by Gear's fifth predictor-corrector method (Haile, 1992) with a time step of t*=0.001 i.e. 2.5 fs . Note that all quantities with an asterisk in this paper, such as L * , R * , ρ * , t*, etc., are non-dimensionalized in terms of σ, , and m, i.e. L*=L/σ, R*=R/σ, ρ * =Nσ 3 /V, t*= t ( /m) 1/2 /σ, T*= k B T/ .  In this research, a cylindrical liquid thread of length L* and radius R* is placed at the center of the computational domain and the remaining space is vacuum. The initial density of the liquid argon is ρ L * =0.819. The system temperature is kept at T*=0.75 or 1.0. These dimensionless values correspond to ρ L =1223 kg/m 3 for argon and T=70K or 93.3K. Note that these two temperatures are below the critical temperature 150K of argon.
The procedure for MD simulation includes three stages : initialization, equilibration and production. Initially, equilibration is performed for liquid argon molecules in a rectangular parallelepiped with length and width equal to the liquid thread diameter (D*=2R*=4, 6, 8) and with height equal to the side length of the computational domain (L*=10, 16, 20, 24, 30, 60, 120, 240, 480). The initial velocities of molecules are decided by normal random numbers. Velocity rescaling is performed at each time step by Eq.
(2) to make sure that the molecules are at the desired temperature T* : where v i new and v i old are the velocities of molecule i after and before correction, respectively, and T D and T A are the desired and the actual molecular temperatures, respectively. The liquid molecules are equilibrated for 10 6 time steps at the desired temperature T*. The achievement of equilibrium state is confirmed by obtaining the radial distribution function.
After the liquid molecules are equilibrated, the rectangular parallelepiped for the liquid molecules is truncated to the desired cylindrical liquid thread by removing unwanted regions. The cylindrical liquid thread then is put into the computational domain and the production stage proceeds. A minimum image method and the Verlet neighbor list scheme (Haile, 1992) to keep track of which molecules are actually interacting at a given time interval of 0.005 are used in the equilibration and the production stages.

Liquid thread vaporization process
In the following discussion, a cylindrical liquid thread of length L* and radius R* is placed at the center of the simulation domain and the remaining space is in vacuum, as illustrated in Fig.1. Simulation domain dimensions, temperatures and number of molecules, are listed in Table 1.
Figures 2~4 show the vaporization processes of liquid threads at T*=0.75 cases 1~23 in Table 1 . Note that L*=120 and R*=3 correspond to L=42.5nm and R=1.06nm, respectively. It is found that when R*=2, the liquid thread breaks up into drops for all fundamental cell lengths. For L*=10, 20, 30 and 60, the liquid thread ruptures only from its two ends, i.e. the top and bottom surfaces of the fundamental cell, and gets shorter due to the contraction motion in its axial direction. Only one liquid particle is formed. On the other hand, for L*=120, 240 and 480, the thread ruptures not only from the top and bottom surfaces of the fundamental cell but also from its interior section. The number of liquid particles produced for L*=120, 240 and 480 is 2, 5 and 10, respectively. During the vaporization process, collision and coalescence of the liquid particles occur. The two liquid particles produced for L*=120 coalesce into one liquid particle. The number of liquid particles at t*=1000 for L*=240 and 480 is 3 and 5, respectively. When R*=4, the liquid thread remains intact for L*=10 and 20. For L*=24, 30, 60, 120 and 240, the liquid thread ruptures merely from its two ends and only one liquid particle is formed. For L*=480, the thread ruptures not only from the top and bottom surfaces of the fundamental cell but also from its interior section. Three liquid particles are produced. The liquid particles eventually coalesce into one liquid particle.  Table 1 . Unlike the liquid threads of the same radius R*=4 at a lower temperature T*=0.75 cases 16~23 in Table 1 , the liquid thread remains intact only for L*=10. For L*=20, 30, 60 and 120, the liquid thread ruptures from its two ends, i.e. the top and bottom surfaces of the fundamental cell, and only one liquid particle is formed. For L*=240 and 480, the thread ruptures not only from the top and bottom surfaces of the fundamental cell but also from its interior section. The number of liquid particles produced for L*=240 and 480 is 2 and 4, respectively. During the vaporization process, collision and coalescence of the liquid particles occur. The liquid particles eventually coalesce into one liquid particle.
www.intechopen.com From the above discussion for Figs.2~5, the following observations can be obtained.
First, if the fundamental cell length is long enough, more than one liquid particle may be produced in the cell. Under such situation, the liquid thread ruptures not only from its two ends but also from the interior of the fundamental cell. On the other hand, if the fundamental cell length is small, the liquid thread may remain intact or produce only one liquid particle in the cell. If the thread breaks up, it ruptures from its two ends only, i.e. the top and bottom surfaces of the fundamental cell, but not from its interior.
Second, a liquid thread with a longer fundamental cell length may produce more liquid particles in the cell.
Third, a thinner liquid thread may produce more liquid particles in the cell.
Fourth, a liquid thread at a higher temperature may produce more liquid particles.

Stability analysis
According to the classical theory by Rayleigh (1879), the radius of the liquid thread, R, is related to the critical wavelength of perturbation, c , as A liquid thread will break up into drops if the axial wavelength of the surface perturbation L＞ c . If L＜ c , the thread is stable and will remain intact. Owing to the periodic boundary conditions in this study, the fundamental cell size L* can be regarded as the longest wavelength of the perturbation.
For R*=2, c * is 12.6. According to Rayleigh's stability criterion, a liquid thread of radius R*=2 and length L*=10 is expected to be stable and remain intact. However, from Fig.21 or Table 1 (Case 1) , it is seen that the liquid thread of R*=2 and L*=10 is unstable, which implies that Rayleigh's stability criterion overpredicts the stable domain as compared to the MD simulation results. Rayleigh's stability criterion holds for the remaining liquid threads with L*＞ c * (L*=20~480) , which are unstable and break up into drop(s).
For R*=3, c * is 18.9. From Fig.3 or Table 1 (Cases 8 and 9) , the liquid thread of R*=3 and L*=10 remains intact while the liquid thread of R*=3 and L*=16 is unstable. Rayleigh's stability criterion holds for the former case but violates the MD simulation results for the latter case. The remaining liquid threads with L*＞ c * (L*=20~480) are unstable and break up into drop(s). Rayleigh's stability criterion holds for these cases.
For R*=4, c * is 25.2. From Fig.4 or Table 1 (Cases 16~23) , the liquid threads of L*=10 and 20 remain intact while the liquid thread of L*=24 is unstable. Rayleigh's stability criterion holds for the first two cases but violates the MD simulation results for the third case. The remaining liquid threads with L*＞ c * (L*=30~480) are unstable and break up into drop(s). Rayleigh's stability criterion holds for these cases.
If the temperature of the liquid threads of R*=4 is increased to T*=1.0, from Fig.5 or Table 1 (Cases 24~30) , only the liquid thread of L*=10 remains intact while, in contrast to the situation of T*=0.75, the liquid thread of L*=20 is unstable. Rayleigh's stability criterion holds for the former case but violates the MD simulation results for the latter case. The remaining liquid threads with L*＞ c * (L*=30~480) are unstable and break up into drop(s). Rayleigh's stability criterion holds for these cases.
Recently, Kim, Lee, Han and Park (2006) proposed a new linear relation from their MD simulation results. From their study, Rayleigh's linear relation can be modified as From Eq.(4), c * for R* of 2, 3 and 4 are 5.9, 13.6 and 21.2, respectively. Then, from Table 1 (Cases 1~23) , it is seen that Kim's stability criterion holds for the liquid threads at T*=0.75. On the other hand, at T*=1.0 (Cases 24~30) , the liquid thread of R*=4 and L*=10 remains intact while the liquid thread of R*=4 and L*=20 is unstable. Kim's stability criterion holds for the former case but violates the MD simulation results for the latter case. The remaining liquid threads with L*＞ c * (L*=30~480) are unstable and break up into drop(s). Kim's stability criterion holds for these cases.
From the above discussion, the following observations can be made.
First, a liquid thread with a longer fundamental cell length is more unstable.
Second, a thinner liquid thread is more unstable.
Third, a liquid thread at a higher temperature is more unstable.
Fourth, the trends of linear stability theories agree with MD simulation results. However, Rayleigh's stability criterion overpredicts stable domain as compared to the MD simulation results. Kim's stability criterion gives more accurate predictions. However, it overpredicts the stable domain at a higher temperature as compared to the MD simulation results.

Density distribution
In practical applications, e.g. combustor or printer, faster vaporization is usually desirable. Criteria have to be made to quantify the discussion regarding the vaporization process of a liquid thread. In this research, a liquid thread i s c o n s i d e r e d t o v a p o r i z e f a s t e r i f t h e distribution of molecules reaches uniform state quicker during the vaporization process. This criterion essentially concerns with the evolution of the density distribution. The density at a specified point in the fundamental cell can be defined as where V is a small volume surrounding the point considered and N is the number of molecules inside the volume V. The density defined by Eq.(5) is actually an averaged density of a small volume surrounding the point considered. The value will approach the density of a specified point if the volume V shrinks to that point. However, for a meaningful density field, the volume V can not be too small because when V becomes too small, it is difficult to obtain a definite value for N/ V. In this study, the volume V is taken to be a sphere with non-dimensionalized radius R*=2 and with its center located at the point considered. This is an optimal choice after numerical test. The time averaged density uniformity factor, f  , in a time interval of t*=0 to 1000, as shown in Fig.6 or listed in Table 1, can be used to indicate the vaporization speed of a liquid thread. The density uniformity factor, f , is defined as where N is the total number of molecules in the fundamental cell, i.e. N mol in Table 1, ρ* and V are the density and volume of molecule i, respectively, as defined by Eq.(5), and ρ* eq is the density value when the molecules are uniformly distributed, i.e. ρ* eq ≡N mol / Vol, where Vol is the volume of the fundamental cell. The density uniformity factor, f , as defined by Eq.(6) represents the deviation from uniform state. From Fig.6 or Table 1, the following observations can be drawn.
First, a liquid thread with a shorter fundamental cell length evaporates quicker. Take the liquid threads of R*=2 (Cases 1~7) as an example. The liquid thread of L*=10 (Case 1) produces one liquid particle while the liquid thread of L*=480 (Case 7) produces ten liquid particles. The length of the latter liquid thread is forty-eight times of the former one. This implies that if forty-eight liquid threads of L*=10 are connected in series, there will be fortyeight liquid particles. It is known that molecular interaction plays an important role in the vaporization process. More liquid particles can provide more molecular interactions and this is conducive to vaporization.
Second, a liquid thread evaporates quicker at a higher temperature.
Third, a liquid thread with a higher f  is more unstable and produces more liquid particles in the fundamental cell. www.intechopen.com

Molecular dynamics analysis of the vaporization process for two nano-scale liquid threads coexisting in a periodic fundamental cell
As mentioned in section 2.1, previous studies of nano-scale liquid threads have mostly been devoted to the investigation of a single liquid thread in a periodic fundamental cell. Because of the interaction between the two nano-scale liquid threads coexisting in a periodic fundamental cell, the vaporization process is different from that of a single liquid thread in a periodic fundamental cell. This section discusses the influences of the liquid thread radius, fundamental cell length, and relative position of the two threads. Snapshots of molecules, the number of liquid particles formed, and density field are analyzed. Two linear stability criteria, namely, Rayleigh's stability criterion and Kim's stability criterion, are accessed for their validity in molecular scale. This approach will be helpful for the understanding and prediction of the atomization process.
In this study, two cylindrical liquid threads of length L* and radius R* are placed in the computational domain and the remaining space is a vacuum. The relative position of the two threads is controlled by L 1 *, L 2 * and L 4 *, as shown in Fig.7. L 3 * is set to be zero, i.e., the two liquid threads are kept at x=0 initially and are shifted only in the y direction. The initial density of the liquid argon is ρ L * =0.819 and the system temperature is kept at T*=0.75. These dimensionless values correspond to ρ L =1223 kg/m 3 for argon and T=70K, which is below the critical temperature (150K) of argon. Simulation domain dimensions, the relative position of the two threads and number of molecules are listed in Table 2, together with some simulation results.  Figure 8 shows the vaporization process of two liquid threads of R*=4, L*=30, L 1 *=10, L 2 *=10 and L 4 *=15 (case 1 in Table 2) . Note that L*=30 and R*=4 correspond to L=10.6nm and R=1.41nm, respectively. It can be seen that the two liquid threads quickly coalesce into a single thread and remain intact. No liquid particles are formed during the vaporization process. Figure 9 shows the vaporization process of two liquid threads of smaller radius R*=2 but with the same length and relative position : L*=30, L 1 *=10, L 2 *=10 and L 4 *=15 (case 2 in Table 2) . For this case in which the liquid threads are thinner, it is seen that the two threads rupture from their two ends, i.e., the top and bottom surfaces of the fundamental cell, and get shorter due to the contraction motion in their axial directions. Two liquid particles are formed and they quickly coalesce into a single particle by collision and coalescence. The coalesced liquid particle prevails during the subsequent vaporization process.

Case
No.
Corresponding  Table 2. Simulation domain dimensions, relative position of the two threads, number of molecules, and simulation results (L 3 * =0) From Table 2, it is seen that the value of f  is smaller for the thinner liquid threads. As stated in the previous section, this implies that the thinner liquid threads evaporate more quickly. In addition, it can be seen from Table 2 that the thinner liquid threads produce more liquid particles. In our previous study (Yeh, 2009b), we found that a single liquid thread is more unstable and produces more liquid particles when it is thinner. The present results for two liquid threads coexisting in a fundamental cell corroborate the findings of the previous study of a single liquid thread.  Table 2) . It is observed that the two threads rupture from their two ends, i.e., the top and bottom surfaces of the fundamental cell, and get shorter due to the contraction motion in their axial directions. Two liquid particles are formed and prevail during the subsequent vaporization process. In contrast to cases 1 and 2 for shorter liquid threads, the two liquid particles do not coalesce. If the two liquid threads get closer to each other, e.g., R*=4, L*=60, L 1 *=10, L 2 *=25 and L 4 *=30 (case 4 in Table 2) , the vaporization process, which is shown in Fig.11, would be very different from case 3. From Fig.11, it is seen that because the two threads are very close to each other, they quickly coalesce into a single thread. The coalesced liquid thread gets shorter due to the contraction motion in its axial direction and finally evolves into a single liquid particle.
From Table 2, it can be seen that the value of f  is larger for case 4. This implies that vaporization is slower when the two liquid threads are close to each other. In addition, Table 2 shows that more liquid particles are formed when the separation of the two threads www.intechopen.com Molecular Dynamics -Studies of Synthetic and Biological Macromolecules 362 is larger (case 3 in Table 2) . In our previous study (Yeh, 2009b), we found that molecular interaction plays an important role in the vaporization process. More liquid particles can provide more molecular interactions, which is conducive to vaporization. Comparison of cases 3, 4 and 5 also reveals this tendency. Case 5, which was investigated in our previous study (Yeh, 2009b), denotes the situation for a single liquid thread of the same radius and length as in cases 3 and 4. From Table 2, it is observed that only one liquid particle is formed for case 5 and the vaporization speed for case 5 is slower than case 3, for which two liquid particles are formed. The present results for two liquid threads coexisting in a fundamental cell corroborate the findings of our previous study for a single liquid thread.  Table 2) . Similar to case 3, the two threads rupture from their two ends and get shorter due to the contraction motion in their axial directions. Two liquid www.intechopen.com Analysis of the Atomization Process by Molecular Dynamics Simulation 363 particles are formed and prevail during the subsequent vaporization process. For two closer liquid threads, e.g., R*=2, L*=60, L 1 *=10, L 2 *=25 and L 4 *=30 (case 7 in Table 2) , the vaporization process, which is shown in Fig.13, is different from case 6. The two threads also rupture from their two ends and get shorter due to the contraction motion in their axial directions. Two liquid particles are formed but they quickly coalesce into a single particle and prevail during the subsequent vaporization process. This vaporization process is similar to but slightly different from that for thicker liquid threads (case 4) . The thinner liquid threads (case 7) evolve into two liquid particles and then quickly coalesce into a single particle while the thicker liquid threads (case 4) quickly coalesce into a single liquid thread and then evolve into a single liquid particle.
From Table 2, it is seen that the value of f  is larger for case 7, i.e., vaporization is slower when the two liquid threads are close to each other. In addition, Table 2 shows that more liquid particles are produced when the separation of the two threads is larger (case 6) . Case 8, which was investigated in our previous study (2009b), denotes the situation for a single liquid thread of the same radius and length as those in cases 6 and 7. From Table 2, it is seen that only one liquid particle is formed for case 8 and the vaporization speed for case 8 is slower than case 6, for which two liquid particles are formed. This corroborates the previous results for thicker liquid threads (cases 3, 4 and 5) . Furthermore, comparing cases 6~8 for thinner liquid threads (R*=2) with cases 3~5 for thicker liquid threads (R*=4) , it can be seen that the value of f  is smaller for thinner liquid threads, which implies that thinner liquid threads evaporate more quickly. This also corroborates the previous results for shorter liquid threads (cases 1 and 2). Figure 14 shows the vaporization process of two liquid threads of R*=4, L*=120, L 1 *=60, L 2 *=30 and L 4 *=60 (case 9 in Table 2) . Similar to case 3 for shorter liquid threads, the two threads rupture from their two ends and get shorter due to the contraction motion in their axial directions. Two liquid particles are formed and prevail during the subsequent vaporization process. For two closer liquid threads, e.g., R*=4, L*=120, L 1 *=40, L 2 *=40 and L 4 *=60 (case 10 in Table 2) , the vaporization process, which is shown in Fig.15, is similar to case 9. Two liquid particles are formed and prevail during the subsequent vaporization process. If the two liquid threads get even closer, e.g., R*=4, L*=120, L 1 *=10, L 2 *=55 and L 4 *=60 (case 11 in Table 2) , the vaporization process, which is shown in Fig.16, is quite different from cases 9 and 10. Similar to case 4 for shorter liquid threads, it is observed that because the two threads in case 11 are very close to each other, they quickly coalesce into a single thread before they separately evolve into liquid particles. The coalesced liquid thread then gets shorter due to the contraction motion in its axial direction and finally evolves into a single liquid particle.
Similar to cases 4 and 7, Table 2 shows that the value of f  is larger for case 11, i.e., vaporization is slower when the two liquid threads are very close to each other. Similarly, it is seen that more liquid particles are produced when the separation of the two threads is larger (cases 9 and 10) . Case 12, which was investigated in our previous study (Yeh, 2009b), denotes the situation for a single liquid thread of the same radius and length as those in cases 9~11. From Table 2, it is seen that only one liquid particle is formed for case 12 and the vaporization speed for case 12 is slower than in cases 9 and 10, for which two liquid particles are formed. This corroborates the previous results for shorter liquid threads (cases 3~8) .  Figure 17 shows the vaporization process of two liquid threads of R*=2, L*=120, L 1 *=60, L 2 *=30 and L 4 *=60 (case 13 in Table 2) . For this case of thinner liquid threads, it is seen that the two threads rupture from their two ends and quickly evolve into ten liquid particles. During the vaporization process, collision and coalescence of the liquid particles occur and finally two liquid particles are formed. For two closer liquid threads, e.g., R*=2, L*=120, L 1 *=40, L 2 *=40 and L 4 *=60 (case 14 in Table 2) , the vaporization process, which is shown in Fig.18, is similar to but slightly different from that of case 13. Nine liquid particles are formed initially and they finally evolve into two liquid particles by collision and coalescence. If the two liquid threads get even closer, e.g., R*=2, L*=120, L 1 *=10, L 2 *=55 and L 4 *=60 (case 15 in Table 2) , the vaporization process, which is shown in Fig.19, is different www.intechopen.com Analysis of the Atomization Process by Molecular Dynamics Simulation 365 from those in cases 13 and 14. Nine liquid particles are formed, but they finally evolve into a single particle by collision and coalescence.  Table 2 shows that the value of f  is larger for case 15, i.e., vaporization is slower when the two liquid threads are very close to each other. Similarly, it can also be seen that more liquid particles are formed when the separation of the two threads is larger (cases 13 and 14) . Case 16, which was investigated in our previous study (Yeh, 2009b), denotes the situation for a single liquid thread of the same radius and length as those in cases 13~15. From Table 2, it is seen that less liquid particles are formed for case 16 and the vaporization speed for case 16 is slower than in cases 13 and 14, in which more liquid particles are formed. This corroborates the previous results (cases 3~12) . Furthermore, comparing cases 13~16 for thinner liquid threads (R*=2) with cases 9~12 for thicker liquid threads (R*=4) , it is seen that the value of f  is smaller for thinner liquid threads, which implies that thinner liquid threads evaporate more quickly. This also corroborates the previous results for shorter liquid threads (cases 1~8). Figure 20 shows the vaporization process of two liquid threads of R*=4, L*=240, L 1 *=120, L 2 *=60 and L 4 *=120 (case 17 in Table 2) . Similar to cases 3 and 9 for shorter liquid threads, the two threads rupture from their two ends and get shorter due to the contraction motion in their axial directions. Two liquid particles are formed and prevail during the subsequent vaporization process. For two closer liquid threads, e.g., R*=4, L*=240, L 1 *=40, L 2 *=100 and L 4 *=120 (case 18 in Table 2) , the vaporization process, which is shown in Fig.21, is similar to www.intechopen.com Molecular Dynamics -Studies of Synthetic and Biological Macromolecules 366 that in case 17. If the two liquid threads get even closer, e.g., R*=4, L*=240, L 1 *=10, L 2 *=115 and L 4 *=120 (case 19 in Table 2) , the vaporization process, which is shown in Fig.22, is quite different from those in cases 17 and 18. Similar to cases 4 and 11 for shorter liquid threads, it is seen that because the two threads in case 19 are very close to each other, they quickly coalesce into a single thread before they separately evolve into liquid particles. The coalesced liquid thread then becomes shorter due to the contraction motion in its axial direction and finally evolves into a single liquid particle. Fig. 19. Vaporization process of two liquid threads of R*=2, L*=120, L 1 *=10, L 2 *=55, L 3 *=0, L 4 *=60 From Table 2, it is observed that the value of f  is larger for case 19, i.e., vaporization is slower when the two liquid threads are very close to each other. Similarly, it can be seen that more liquid particles are produced when the separation of the two threads is larger (cases 17 and 18) . Case 20, which was investigated in our previous study (Yeh, 2009b), denotes the situation for a single liquid thread of the same radius and length as those in cases 17~19. From Table 2, it is seen that only one liquid particle is formed for case 20 and the vaporization speed for case 20 is slower than those in cases 17 and 18, for which two liquid particles are formed. This corroborates the previous results for shorter liquid threads (cases 3~16) . If the length of the liquid thread is further increased, e.g. R*=4 and L*=480, more liquid particles are formed. Figure 23 shows the vaporization process of two liquid threads of R*=4, L*=480, L 1 *=160, L 2 *=160 and L 4 *=240 (case 21 in Table 2) . Eight liquid particles are formed and they evolve into four liquid particles by collision and coalescence. For two closer liquid threads, e.g., R*=4, L*=480, L 1 *=40, L 2 *=220 and L 4 *=240 (case 22 in Table 2) , the vaporization process, which is shown in Fig.24, is similar to but slightly different from that in case 21. Six liquid particles are formed initially and they finally evolve into two liquid particles by collision and coalescence. If the two liquid threads get even closer, e.g., R*=4, L*=480, L 1 *=10, L 2 *=235 and L 4 *=240 (case 23 in Table 2) , the vaporization process, which is shown in Fig.25, is very different from those in cases 21 and 22. Because the two threads in case 23 are very close to each other, they quickly coalesce into a single thread before they separately evolve into liquid particles. The coalesced liquid thread ruptures from its ends and interior and gets shorter due to the contraction motion in its axial direction and finally evolves into a single liquid particle by collision and coalescence. The vaporization process for case 23 (Fig.25) is similar to but slightly different from that in case 19 (Fig.22) . For both cases, the two threads quickly coalesce into a single thread because they are very close to each other. However, the coalesced liquid thread in case 19 evolves into a single liquid particle without rupturing in its interior while the coalesced liquid thread for case 23 ruptures in its interior before it finally evolves into a single liquid particle.  Table 2, it is observed that the value of f  is larger for case 23, i.e., vaporization is slower when the two liquid threads are very close to each other. Similarly, it can be seen that more liquid particles are produced when the separation of the two threads is larger (cases 21 and 22) . Case 24, which was investigated in our previous study (Yeh, 2009b), denotes the situation for a single liquid thread of the same radius and length as those in cases 21~23. From Table 2, it can be seen that less liquid particles are formed for case 24 and the 368 vaporization speed for case 24 is slower than those in cases 21 and 22, in which more liquid particles are formed. This corroborates the previous results for shorter liquid threads (cases 3~20). From the above discussion, the following observations can be made.
First, if the fundamental cell length is small, the liquid threads may remain intact or evolve into only one liquid particle in the cell. If the threads break up in this case, they rupture from their ends only, i.e. the top and bottom surfaces of the fundamental cell, but not from their interiors. On the other hand, if the fundamental cell length is larger, more than one liquid particle may be produced in the cell and the liquid threads may rupture not only from their ends but also from their interiors.
Second, thinner liquid threads may produce more liquid particles in the cell and evaporate more quickly.
Third, more liquid particles are formed when the separation of the two threads is larger. Moreover, vaporization is slower when the two liquid threads are close to each other.
Fourth, on the basis of identical liquid thread radius and length, liquid threads that produce more liquid particles evaporate more quickly.

Stability analysis
In the following discussion, two linear stability criteria, namely Rayleigh's stability criterion (1879) and Kim's stability criterion (2006), are accessed for their validity in molecular scale. This is the first study to investigate the instability of two nano-scale liquid threads coexisting in a periodic fundamental cell by MD simulation.
According to classical theory by Rayleigh (1879), the radius of a liquid thread, R, is related to the critical wavelength of perturbation, c , by Eq.(3) in section 2.1.3. A liquid thread will break up into drops if the axial wavelength of the surface perturbation L＞ c . If L＜ c , the thread is stable and will remain intact. Owing to the periodic boundary conditions in this study, the fundamental cell size L* can be regarded as the longest wavelength of the perturbation.
For R*=4, c * is 25.2. According to Rayleigh's stability criterion, liquid threads of radius R*=4 and length L*=30, 60, 120, 240 and 480 are expected to be unstable and will break up into drop(s). However, from Figs.8, 10, 11, 14~16, 20~25 or Table 2 (cases 1, 3~5, 9~12, 17~24) , the liquid threads of R*=4 and L*=30 (case 1) remain intact while the remaining liquid threads are unstable and break up into drop(s). Rayleigh's stability criterion violates the MD simulation result for case 1 but holds for the remaining cases. The invalidity of Rayleigh's stability criterion for case 1 can be interpreted according to the vaporization process shown in Fig.8. It is seen that the two liquid threads quickly coalesce into a single thread and remain intact. This is somewhat like combining the two threads to form a thread of a larger radius, which can increase the critical wavelength of perturbation, c , according to Eq.(3) and, therefore, broaden the stable domain.
Recently, Kim, Lee, Han and Park (2006) proposed a new linear relation from their MD simulation results. From their study, Rayleigh's linear relation can be modified by Eq.(4) in section 2.1.3. From Eq.(4), c * for R* of 2 and 4 are 5.9 and 21.2, respectively. Then, from Table  2, it is seen that Kim's stability criterion violates the MD simulation result for case 1 while it holds for the remaining cases. The validity of Kim's stability criterion is similar to that of Rayleigh's stability criterion.
From the above discussion, it is found that the trends of Rayleigh's stability criterion and Kim's stability criterion agree with MD simulation results. However, when the two liquid threads coalesce into a single thread and remain intact, the critical wavelength of perturbation may be increased and the stable domain is broadened. In such a situation, Rayleigh's stability criterion and Kim's stability criterion underpredict the stable domain.

Molecular dynamics simulation for the atomization process of a nanojet
In this section, liquid argon nanojets made of 44000 Lennard-Jones molecules are investigated under various simulation parameters to examine their influence on the nanojet atomization process. Snapshots of the molecules, evolution of the density field, and evolution of the intermolecular force are analyzed. This can provide insight into the fundamental mechanism of the atomization process and will be helpful for the design of nanojet devices such as nano-printer or nano-sprayer.
The nano-atomizer is schematically shown in Fig.26. The simulation domain comprises a cubical box of side length 3600, with periodic boundary conditions applied in all three directions. The nano-atomizer is placed at the center of the box. Simulation parameters are listed in Table 3, which include nano-atomizer dimensions, temperatures, number of molecules and simulation results.
The argon molecules inside the nano-atomizer are liquid and the nano-atomizer is made of rigid argon molecules. A push panel composed of 600 argon molecules is constructed with a downward velocity of 120m/s. In this study, the interactions among liquid argon molecules, nano-atomizer and push panel are taken into account.   Table 3. Nano-atomizer dimensions, temperatures, number of molecules and simulation results In the following discussion, a liquid argon nanojet of length L 5 * and diameter 2L 3 * +D * is pushed by a panel into vacuum through a nano-nozzle of orifice diameter D * , as illustrated in Fig.26. Simulation conditions are listed in Table 3. Figure 27 shows the atomization process for a nanojet of L 2 * =L 3 * =L 4 * =5.73, L 5 * =76.7, D * /2=8.81 and T*=0.75, which corresponds to a nanojet of length 26.2 nm and diameter 10 nm, and a nano-nozzle of orifice length 2nm, diameter 6 nm, as well as an actual temperature of 70 K. The dot in Fig.27 indicates the center of the molecule. From the figure it is found that the nanojet does not break up. Owing to the low temperature, the molecular kinetic energies are so low that the molecules congregate near the orifice exit. Very few liquid molecules evaporate at this low temperature. At a higher temperature T*=1.5 (140 K) , as depicted in Fig.28, the molecules leave the orifice exit earlier than at T*=0.75, due to their higher molecular kinetic energies. More liquid molecules evaporate at this higher temperature. However, like at T*=0.75, the nanojet has not broken up before t*=80. If the temperature is further increased to T*=2.0 (187 K) , as shown in Fig.29, evident evaporation is observed. Many evaporated molecules are produced and the non-evaporated liquid molecules concentrate within the central region. Figure 30 depicts the snapshots for temperature T*=3.0 (278 K) . It is observed that breakup of the nanojet occurs and its spray angle is larger than at T*=2.0. The spurted molecules from the nano-atomizer are more evenly distributed at this temperature. If the temperature is further increased to T*=4.5 (420 K) , as shown in Fig.31, the spray angle is even larger than at T*=3.0 and the spurted molecules from the nano-atomizer are much more uniformly distributed as compared to the lower temperature cases. Comparison of Figs.27,28,29,30 and 31 reveals that the liquid nanojet evaporates quicker at higher temperatures. This will be further illustrated in later sections discussing the density distribution and the intermolecular force.

Nanojet atomization process
www.intechopen.com  Figure 52 reveals similar tendency. Note that the nano-atomizer in this research is basically a plain-orifice atomizer. As pointed out by Lefebvre (1989) , in a practical plain-orifice atomizer, resistance increases with nozzle orifice length/diameter ratio. Therefore, a nanojet from a nozzle with a smaller orifice length/diameter ratio moves farther due to its smaller resistance. This will be further illustrated in later sections discussing the density distribution and the intermolecular force.

Density distribution
It is important that the system be in equilibrium state before statistical values of the local properties can be taken. However, owing to the computational capacity limitations, the MD simulation can not proceed to a macroscopically long period. Nevertheless, the purpose of this paper is not to discuss statistical values of the local properties but to investigate the atomization process of a nanojet, which is important and conducive to the understanding of the fundamental mechanism of the atomization process. Criteria have to be made to quantify the discussion regarding the nanojet atomization process. Unfortunately, such criteria are still arbitrary in the literature. Because the system temperature in this study is kept at the desired temperature, a constant temperature criterion is not suitable for the discussion of the atomization process. In this research, a nanojet is considered to vaporize faster if the distribution of molecules reaches a uniform state quicker during the atomization process. This criterion essentially concerns with the evolution of the density distribution. The density at a specified point in the fundamental cell can be defined as Eq.(5). In this study, the volume δV is taken to be a sphere with nondimensionalized radius R*=2 and with its center located at the point considered. This is an optimal choice after numerical test.  Figure 34 shows the evolution of density uniformity factor for nanojets at different temperatures and the conditions of L 2 * =L 3 * =L 4 * =5.73, L 5 * =76.7, D * /2=8.81 ( cases 1~5 in Table 3 ). The density uniformity factor is defined as Eq.(6). From Fig.34 it is observed that a higher temperature nanojet evaporates faster than a lower temperature one and this corroborates the results of Figs.27~31 as discussed in section 2.3.1. The time averaged value of the density uniformity factor, f  , in a time interval of t*=0 to 80, as listed in Table 3, also reveals this observation. In Fig.34, it is noted that at lower temperatures (T*=0.75 and 1.5) , the density uniformity factor increases first and then decreases. For a lower temperature nanojet, the momenta of the liquid molecules away from the push panel in the nano-atomizer are low while the molecules near the push panel have relatively higher momenta due to the action of the push panel. This results in a compression effect that leads to the increase of the density uniformity factor at the earlier stage of the atomization process; while at a later stage, the density uniformity factor drops because of the ejection of the molecules.  Table 3) Figure 35 shows the evolution of density uniformity factor for different orifice lengths on the basis of identical nozzle orifice diameter (6 nm) and the conditions of L 3 * =5.73, L 5 * =76.7, T*=2.0 ( cases 3, 6 and 7 in Table 3 ). It is observed that a nanojet with a shorter orifice length evaporates quicker. The time averaged value of the density uniformity factor, f  , in a time interval of t*=0 to 80, as listed in Table 3, also reveals this observation. This corroborates the results of Figs.32(a), (b) and (c) as discussed in section 2.3.1. In addition, it is also observed that, as time elapsed, the influence of the orifice length mitigates. This is because as time elapsed, more and more molecules spurt from the atomizer and hence the interaction between the liquid molecules and the rigid atomizer molecules mitigates due to the decrease of number of molecules inside the atomizer. Figure 36 also reveals this tendency. However, as can be observed from Fig.36, the influence of the orifice length becomes less pronounced at a higher temperature because of the higher molecular kinetic energy to overcome the resistance caused by the orifice. Figure 37 shows the evolution of density uniformity factor for different orifice diameters on the basis of identical nozzle orifice length (2 nm) and the conditions of L 4 * =5.73, L 5 * =76.7, T*=2.0 (cases 3 and 8 in Table  3). It is observed that a nanojet with a larger orifice diameter evaporates quicker. The time averaged value of the density uniformity factor, f  , in a time interval of t*=0 to 80, as listed in Table 3, also reveals this observation. This corroborates the results of Figs.32(a) and (d) and also reveals previous observation that a nanojet from a nozzle with a smaller orifice length/diameter ratio evaporates quicker. Similar tendency is obtained from Fig.38 for a higher temperature nanojet.

Intermolecular force
The intermolecular force is an indication of the surface tension experienced by the liquid particles and has a great effect upon the atomization process. Lefebvre (1989), Chigier (1999) and Hiroyasu (2000) pointed out that surface tension and interfacial force are the major controlling mechanisms for atomization. Owing to the vacuum environment, the aerodynamic effect on the atomization process is negligible in this study. Thus, the surface tension becomes the major controlling mechanism for the atomization process. Figure 39 shows the evolution of averaged non-dimensionalized intermolecular force for nanojets with different temperatures and the conditions of L 2 * =L 3 * =L 4 * =5.73, L 5 * =76.7, D * /2=8.81 ( cases 1~5 in Table 3 ). The averaged non-dimensionalized intermolecular force at time t* is defined as where N is the total number of molecules in the fundamental cell and F * i,t* is the resultant force of the non-dimensionalized intermolecular force vector acting on molecule i at time t* , i.e. F * i,t* = ( F * x,i,t* 2 + F * y,i,t* 2 + F * z,i,t* 2 ) 1/2 , where F * x,i,t* , F * y,i,t* and F * z,i,t* are the components of the intermolecular force vector at the x, y and z directions, respectively, acting on molecule i at time t*. Note that in the above definition of * t* F , N is the total number of molecules in the fundamental cell, which includes liquid, vapor and solid molecules (atomizer and push panel) ; while in the definition of density uniformity factor, Eq.(6), N is only the initial number of liquid molecules in the fundamental cell, i.e. the solid molecules are excluded. The intermolecular force diminishes with time because of the increase of distances between molecules as the nanojet vaporizes. From Fig.39, it is observed that a higher temperature nanojet evaporates faster than a lower temperature one. This corroborates the results of Figs.27~31 discussed in section 2.3.1 and Fig.34 in section 2.3.2. In Fig.39, it is also noted that although a higher temperature nanojet has a larger intermolecular force at the earlier stage of the atomization process due to its higher momentum, it evaporates faster and therefore the intermolecular force decays quicker. Figure 40 shows the time averaged value of the averaged non-dimensionalized intermolecular force for different orifice lengths on the basis of identical nozzle orifice diameter (6 nm) and the conditions of L3 * =5.73, L 5 * =76.7, T*=2.0 ( cases 3, 6 and 7 in Table 3). It is observed that a nanojet with a shorter orifice length evaporates quicker. This corroborates the results of Figs.32(a), (b) and (c) as discussed in section 2.3.1 and Fig.35 as discussed in section 2.3.2. In addition, it is also observed that, as time elapsed, the influence of the orifice length mitigates. As explained in section 2.3.2, more and more molecules spurt from the atomizer as time elapsed. This causes the interaction between the liquid molecules and the rigid atomizer molecules to mitigate due to the decrease of number of molecules inside the atomizer. Figure 41 also reveals this tendency. However, as can be observed from Fig.41, the influence of the orifice length becomes less pronounced at a higher temperature because of the higher molecular kinetic energy to overcome the resistance caused by the orifice. Figure 42 shows the evolution of averaged non-dimensionalized intermolecular force for different orifice diameters on the basis of identical nozzle orifice length (2 nm) and the conditions of L 4 * =5.73, L 5 * =76.7, T*=2.0 ( cases 3 and 8 in Table 3 ). It is observed that a nanojet with a larger orifice diameter evaporates quicker. This corroborates the results of Figs.51(a) and (d) as discussed in section 2.3.1 and Fig.37 as discussed in section 2.3.2 and also reveals that a nanojet from a nozzle with a smaller orifice length/diameter ratio evaporates quicker. Similar tendency can be observed from Fig.43 for a higher temperature nanojet.  Table 3 )

Conclusions
The purpose of this study is to provide a thorough investigation and understanding of the microscopic evolution and detailed mechanism of the flow inside and outside the atomizers. The main findings of this study are listed below.
1. For a single nano-scale liquid thread in a periodic fundamental cell, the liquid thread is more unstable and produces more liquid particles in the fundamental cell when it is thinner or at a higher temperature. In addition, a liquid thread with a longer fundamental cell length is also more unstable and produces more liquid particles in the fundamental cell, but it evaporates slower. The trends of linear stability theories agree with MD simulation results. However, Rayleigh's stability criterion overpredicts stable domain as compared to the MD simulation results. Kim's stability criterion gives more accurate predictions but overpredicts the stable domain at a higher temperature. Finally, a liquid thread with a higher time averaged density uniformity factor, f  , is more unstable and produces more liquid particles in the fundamental cell. 2. For two nano-scale liquid threads coexisting in a periodic fundamental cell, the two liquid threads may remain intact or evolve into only one liquid particle if the fundamental cell length is small. If the threads break up in this case, they rupture from their ends only, i.e., the top and bottom surfaces of the fundamental cell, but not from their interiors. On the other hand, if the fundamental cell length is larger, more than one liquid particle may be produced in the cell and the liquid threads rupture not only from their ends but also from their interiors. It is also found that thinner liquid threads may produce more liquid particles in the cell and evaporate more quickly. In addition, more liquid particles are formed when the separation of the two threads is larger. Moreover, vaporization is slower when the two liquid threads are close to each other. On the basis of identical liquid thread radius and length, liquid threads that produce more liquid particles evaporate more quickly. Finally, the trends of Rayleigh's stability criterion and Kim's stability criterion agree with MD simulation results. However, when the two threads coalesce into a single one and remain intact, the critical wavelength of perturbation may be increased and the stable domain is broadened. In such a situation, Rayleigh's stability criterion and Kim's stability criterion underpredict the stable domain. 3. For the atomization process of a nanojet, it is found that a liquid nanojet evaporates faster at a higher temperature. On the basis of identical nozzle orifice diameter, a nanojet from a nozzle with a shorter orifice length evaporates quicker. However, the influence of the orifice length mitigates as time elapsed. In addition, the influence of the orifice length becomes less pronounced at a higher temperature. On the other hand, on the basis of identical nozzle orifice length, a nanojet from a nozzle with a larger orifice diameter evaporates quicker. The present simulation results reveal that a nozzle with a smaller orifice length/diameter ratio produces better atomization. This corroborates the results from conventional macroscopic analysis.

Acknowledgment
The author gratefully acknowledges the grant support from the National Science Council, Taiwan, R.O.C., under the contract NSC100-2221-E-150-047.