Jtam-A4.dvi


JOURNAL OF THEORETICAL

AND APPLIED MECHANICS

54, 3, pp. 881-892, Warsaw 2016
DOI: 10.15632/jtam-pl.54.3.881

ANALYSIS OF PRESSURE BEHAVIOR IN A TEMPERATURE

CONTROLLED MOLECULAR DYNAMIC FLOW

Hamed R. Najafi, S.M. Hossein Karimian

Department of Aerospace Engineering, Amirkabir University of Technology, Tehran, Iran

e-mail: hamedramezani@aut.ac.ir; hkarim@aut.ac.ir

Thermo-fluid properties are required for numerical modeling of nano/micro devices. These
properties aremostly obtained fromresults ofmolecular dynamics (MD) simulations.There-
fore, efforts havebeenput indevelopingmethods for numerical evaluationof fluidproperties,
such as pressure. In this paper, the pressure behavior in a controllable nanochannel flow is
investigated.The nanoflowfield is created by imposing a gradient of amacroscopic property
such as temperature. Details of the pressure calculationmethod in a molecular system and
its sensitivity to the approximationsmade are described first. The effect of temperature rise
in a uniform flow on the pressure field is studied next. Then, in the flowunder a fixedmean
velocity condition, the effect of temperature gradient as a controllable property on the pres-
sure field of nanoflow is studied. Velocity, pressure andmolecular density of nanoflowswith
various temperature gradients and different temperature levels are investigated as well. It
has been found that the temperature level at which the temperature gradient is imposed, is
important. A fixed temperature gradient will not always lead to the same pressure gradient
at different temperature levels. Furthermore, quite interestingly, it is observed that at a fixed
temperature gradient, with the variation of mean velocity the pressure field also varies.

Keywords: nanofluid, molecular dynamics, pressure, bin size, sampling, periodic flow

1. Introduction

Flow behavior at micro and nano scales has been a subject of interest in the recent years. The
flow in a nanochannel, as a typical reduced-size fluid flow system, can demonstrate different
aspects of nanofluid systems. This is why many studies can be found in literature focusing on
the simulation of nanochannel flows (Mi and Chwang, 2003).
Molecular Dynamics (MD) is a deterministic method to calculate position of molecules and

their dynamicproperties.Differentpotential functionshavebeen introduced to representmolecu-
lar forces inMDsimulations (see for example Leach, 1999; Sadus, 2002). Extractingmacroscopic
properties such as velocity, temperature andpressure frommicroscopic datahas beenalso a chal-
lenging issue for scientists (Allen and Tildesley, 1987; Karniadakis et al., 2002; Karimian et al.,
2011). In earlier studies, the boundary of the domainwas not of primary importance and simple.
Periodic boundary conditions were used (Stillinger andRahman, 1974; Travis and Evans, 1977;
Koplik et al., 1989; Travis and Gubbins, 2000; Nagayama and Cheng, 2004). Implementation
of boundary conditions in MD simulation has been the center of attention in the recent years,
especially in the field of mechanics (Sun and Ebner, 1992; Huang et al., 2004, 2006; Hanasaki
and Nakatani, 2006). In addition to the inlet/outlet boundary conditions, different approaches
exist in the literature to implement wall boundary conditions inmolecular dynamics. Themost
usual choice is to model walls by two rows of molecules. Wall molecules are allowed to oscillate
about their initial positions at which they are fixed. This is done either by assigning heavy
weights to the wall molecules (Koplik et al., 1988), rescaling the velocity of wall molecules to
their initial values (Huang et al., 2004) or using fictitious springs (Fan et al., 2002; Huang et



882 H.R. Najafi, S.M.H. Karimian

al., 2006; Sofos et al., 2009; Branam andMicci, 2009; Kamali andKharazmi, 2011). In addition
to these, some studies can be found where the walls were modeled differently. For example, see
thework of Ziarani andMohamad (2005), where thewalls weremodeled as reflective. Sincewall
boundary conditions are not the subject of this paper, interested scientists are encouraged to
read the above-mentioned references.

Boundary conditions at the inlet and outlet have beenalso a subject of research.With proper
boundary conditions, the desired flow field can be created within a nanochannel. There are
differentmethods in the literature to create controllable nanoflow. Thesemethods are classified
in four categories. Creating a flowwith an external force or acceleration tomovemolecules in a
specified direction is the most simple and straightforward method. Note that the extra energy
added to the flowfield due to external forcemust be taken out of the domain to hold the energy
balance (Koplik et al., 1988; Fan et al., 2002; Ziarani and Mohamad, 2005; Sofos et al., 2009;
Kamali and Kharazmi, 2011).

Another method found in the literature makes use of motion of a piston or plate upstream
of the inlet to create the flow. In this method, the solution domain is extended to provide extra
space between the piston and the inlet. To control flow condition in the outlet aswell, sometimes
a piston and the required extra space is added to the end of the solution domain. In this case,
there is no need to implement periodic boundary conditions or insert and delete molecules at
the inlet or outlet of the solution domain. This approach works well for dense systems. Since
additional molecules are to be taken into account in these added spaces, extra computational
effort is required in these methods. The other disadvantage is that the solution is limited to
the time in which the pistonmoves in the space between its initial position and the inlet/outlet
(Huang et al., 2004, 2006; Hanasaki and Nakatani, 2006).

The thirdmethod involves using a temperature gradient to create a nanoflow. In this appro-
ach, hot and/or cold walls are embeddedwithin the domain along the flowdirection (Han, 2008;
Liu and Li, 2010; Darbandi et al., 2011). The role of the temperature gradient is to change the
kinetic energy of the fluid. At high temperatures, the kinetic energy and therefore the pressure
of the fluid increase. This pressure gradient causes the fluid to move downstream. The energy
balance can be established between hot and cold plates to prevent accumulation of extra energy
in the domain. It should benoted that in practice, hot or cold temperaturewalls cannot be easily
mountedwithin the flow and, therefore, implementation of this approach is not straightforward.

The fourthmethod creates a flow bywall motion. In this approach, walls of a channel move
in opposite or parallel directions with a constant velocity to create the flow (Zhang et al., 2009;
Kim et al., 2010).

While all of the above-mentioned methods achieve the goal of creating a molecular flow in
nanochannels, they fail to offer a proper approach to produce the predefined flow. The goal
of the present research is to investigate pressure behavior in a controllable nanochannel flow
created by a gradient of a macroscopic property such as temperature. In the following sections,
we first introduce the calculationmethod used here to evaluate pressure in a nanoflow. Thenwe
will discuss details about the nanoflows controlled by the temperature gradient, especially the
behavior of pressure.

2. Molecular dynamics

In molecular dynamics, positions of molecules are determined using Newton’s second law. In-
termolecular fluid-fluid forces are calculated using the Lennard-Jones 12-6 potential equation
(Rapaport, 2004)

Uij =4ε
[( σ

rij

)12
−

( σ

rij

)6]

(2.1)



Analysis of pressure behavior in a temperature controlled... 883

In this equation,Uij is the Lennard-Jones potential, rij is the distance between two interacting
molecules i and j, ε is the energy scale, and σ is the finite distance at which the intermole-
cular potential is zero. The derivative of the L-J potential with respect to rij represents the
intermolecular force

fij =−
dUij

drij
fij =

48ε

r2ij

[( σ

rij

)12
−

1

2

( σ

rij

)6]

rij (2.2)

In the above equation, fij is the force of molecule j acting on molecule i. It is common to
enforce this equation within a cut-off distance to reduce the computational cost. We chose a
cut-off distance of rc = 2.5σ herein, which is usually used in other references to simulate the
argon flow (Koplik et al., 1989; Priezjev, 2007; Sofos et al., 2009).

3. Pressure formula

In amacroscopic analysis, pressure is calculated solely based on the virial equation of state. In a
microscopic systemhowever, pressure is calculated from the virial equation of state on all atoms
plus summation of intermolecular forces multiplied by corresponding distances between them.
For pair-additive force fields, the fluid pressureP can be estimated through the virial equation
of state given below (Allen and Tildesly, 1987)

P =
NkBT

V
+
1

3V

N
∑

i=1

N
∑

j>i

rijfij (3.1)

whereN is the number of molecules in the computational domain, kB is Boltzmann’s constant,
and V and T are volume and temperature of the computational domain, respectively. The first
term in Eq. (3.1) is called the kinetic term and contains temperature of the computational
domain. In a microscopic view, temperature is defined as

T =
m

3NkB

N
∑

n=1

3
∑

i=1

(Vni −V i)
2 (3.2)

wherendenotes themoleculenumber in thedomain, i=1,2,3denote thex,y, andz components
of the atomic velocity Vn, V i is the i-th component of the mean flow velocity V, andm is the
atom mass. The mean velocity can be obtained from any of the several averaging methods like
SAM (Tysanner and Garcia, 2004), CAM (Tysanner and Garcia, 2005) or SMC (Karimian et
al., 2011) that reduces statistical errors in its calculation. In this paper, the SMC method is
used to calculate the mean velocity. As can be observed in Eq. (3.2), temperature is related to
the kinetic energy of molecules in the microscopic view, which is always positive. The velocity
difference within parentheses that represents instantaneous velocity due to thermal fluctuations
is called the virial velocity.
The second term in Eq. (3.1) includes the effect of forces between the molecules; this term

is called the potential term. The kinetic term is always positive but the potential term may be
positive or negative. A zero force exists between twomolecules at a distance of 21/6 based on the
L-J potential model. If the distance between two molecules is higher than this value, the sign
of the pair-additive force will be negative and the molecules attract each other. On the other
hand, when the distance between two molecules is less than this value, the force is of repellant
nature.
The pressure calculation method defined by Eqs. (3.1) and (3.2) is implemented in an in-

house, MD code to determine pressure values from the results of MD simulations. Please note
that this code is used in authors’ previous works in which its accuracy in calculation of ma-
croscopic properties were demonstrated (Karimian et al., 2011; Karimian and Namvar, 2012;
Namvar and Karimian, 2012; Karimian and Izadi, 2013).



884 H.R. Najafi, S.M.H. Karimian

4. Pressure calculation

The purpose of this Section is to validate the pressure calculation. According to Eq. (3.1),
pressure has two parts, kinetic and potential terms. Our experience shows that precautions
shouldbe taken in calculation of pressure inMDsimulation.This is demonstratedby considering
a cubic regionwith itsmolecules in equilibrium,which is a constantNumber,VolumeandEnergy
(NVE) simulation. We know that pressure would be the same all over this region. The size of
this cube is L = 45Å and it contains N = 1728 mono atomic noble gas molecules of Argon.
Atomic diameter of Argon is σ=3.4Å and its energy parameter is ε=1.67·10−21J. A periodic
boundary condition is applied in all directions.Argonmolecules are initially arranged in a lattice
form of Face-Centered Cube (FCC). MD simulation of Argon molecules starts from the initial
temperature of 120K, and the mean velocity of zero. The equation of motion is integrated
over time with ∆t= 10−15 s, using Verlet’s scheme (Verlet, 1967). The number of molecules is
constant during the simulation.Macroscopic properties at any point of themolecular domain is
extracted via sampling and averaging of molecular properties within the control region around
that point, called a bin (Kamiadakis et al., 2002). Four bins with different sizes are selected
for sampling and averaging of flow properties. As shown in Fig. 1, the largest bin covers the
whole cube and contains all argon molecules. Other bins are smaller cubes with lateral sizes of
3/4, 1/2 and 1/4 of the lateral size of the whole domain. In volumetric fraction, these cubes
are 27/64, 1/8, and 1/64 of the whole domain, respectively. Note that the periodic boundary
condition is applied on the cube boundaries. Therefore, similar to other three bins, bin 1/1 is
also simulated as if it is in the middle of the domain. Obviously, pressure obtained from all bin
sizes should be equal. The molecular dynamics simulation is continued for 30000 time steps of
10−15 s, when equilibrium is reached. Although pressure as a macroscopic property should be
independent of the bin size, the results show that different pressure values are obtained. The
origin of this difference is in the potential part of the pressure equation.

Fig. 1. Bin sizes for pressure calculation in a periodic stationarymolecular flow

The kinetic part of pressure, which is solely based on the properties of atoms and is indepen-
dent of molecular interactions, converges to the value of 27MPa in all bins with different sizes.
The potential part of pressure, however, converges to different values of −22.3, −18.4, −1.2,
and 5.9MPa in bins with sizes of 1/1, 3/4, 1/2, and 1/4, respectively. This part of pressure is
calculated from summation of all pair-additive forces between molecules multiplied by the di-
stance between them. In this approach, no cut-off distance is considered for themolecules.More
numerical tests have shown that inclusion of the cut-off distance cannot resolve this problem of
pressure calculation, unless the influence of themolecules outside the calculation bin within the
cut-off distance from its boundaries is taken into account. Having implemented this, calculated
pressures in all of the bins converge to the value of 31.4MPa at the equilibrium state withminor
pressure oscillations and differences in bins with sizes of 1/4 and 1/2 due to the bin size. Note



Analysis of pressure behavior in a temperature controlled... 885

that the number of molecules within a bin plays an important role in calculation of pressure.
More details about the effect of bin size on sampling and averaging can be found in (Karimian
and Izadi, 2013). Calculation of pressure in bin 1/1 with LAMMPS (Plimpton, 1995) converged
to the same value with less than a 3 percent difference, which verifies the present method for
calculation of pressure inMD simulation.

5. Control of pressure by temperature

As noted, macroscopic properties can be extracted from microscopic properties of molecules,
calculated by MD simulation. In contrast to this, here we would like to set a macroscopic
property in the solution domain. For this purpose, consider a periodic flow in a nanochanel
with length of L = 144.6Å and cross-section of 36.15Å×36.15Å. This flow contains N = 940
mono-atomic noble gas molecules of Argon. Periodic boundary conditions are applied on the
boundaries of the domain. Argon molecules are initially arranged in a FCC lattice form, and
solution starts from the initial temperature and velocity of 120K and 5m/s, respectively. As
shown in Fig. 2, the solution domain is divided into 8 equal bins.We are going to set a desired
pressure within the domain and study its effects on the solution; but before that, the solution
is proceeded for 400000 time steps of ∆t = 10−15 s to reach equilibrium. During the solution,
temperature of 120K and mean velocity of 5m/s are continuously enforced in the first bin.
At each time step, the mean velocity in this bin is updated by adding the difference between
desired and calculatedmean velocities to the velocity of eachmolecule inside thebin. In a similar
fashion, at each time step, temperature in the bin is updated by rescaling the virial velocity of
molecules. As expected, temperature and flow velocity within the whole solution domain reach
120K and 5m/s respectively, at the equilibrium. The averaging method of SMC (Karimian et
al., 2011) is applied to calculate pressure usingEq. (3.1) at each time step. Solution convergence
for velocity and pressure in all 8 bins is shown in Fig. 3.

Fig. 2. Periodic flow including 940molecules with 8 bins along the x direction

As expected, the mean velocity has converged to a value of 5m/s in all of the bins. The
convergence of themean velocity values is faster than that for the pressure.We believe that this
is because themean velocity is weakly dependent on interactions betweenmolecules. In Fig. 3b,
pressure has converged to a value of 2.9MPa.
Nowwewould like to set a newvalue for pressure in the solution domain. For a fixed number

ofmolecules, pressure can be set to a newvalue by changing temperature ofmolecules. As noted
above, temperature can be set by rescaling the virial velocity of molecules. To examine this,
temperature is raised up from 120K to 150K in the first bin of this periodic flow at time
4 ·10−10 s; i.e. in the time step of 400001. Then solution is continued for an extra 400000 time
steps. The enforced mean velocity is still 5m/s in the first bin. Analysis of the solution shows
that pressure increases in the first bin and consequently propagates to the rest of the periodic



886 H.R. Najafi, S.M.H. Karimian

Fig. 3. Solution convergence in bins for: (a) velocity in a periodic flowwith enforced velocity of 5m/s
in bin 1, and (b) pressure in a periodic flowwith enforced temperature of 120K in bin 1

Fig. 4. Convergence of pressure in different bins of the periodic flowwith the mean flow velocity
of 5m/s. A raise in temperature from 120K to 150K is enforced in the first bin at the time step 400001

flow until the equilibrium is reached. The convergence of pressure in all of the bins is shown in
Fig. 4. At the equilibrium, pressure reaches 6.1MPa all over the solution domain. Again, the
mean flow velocity also reaches the value of 5m/s.

In thenext step,wewould like to seehoweffectively the temperaturegradient alongaperiodic
flow can alter pressure, or how quick would be the response of pressure to the temperature
gradient along the flow. In the first 400000 time steps, the flow field reaches its equilibrium for
the enforced temperature and mean flow velocity of 120K and 5m/s, respectively, in the first
bin. Then, at the time step 400001, a temperature difference of 30K is applied between bins 3
and 7, by setting temperature values to 180K and 150K, respectively.

During the first 400000 time steps, SMC averaging uses all solution data from the first time
step. However, after the implementation of temperature change at the time step 400001, we
cannot includemuch data from the previous time steps in SMCaveraging. If we start SMCave-
raging promptly from the time step 400001, undesirable pressure fluctuations will appear in the
results. In this case, noticeable extra time steps would be required to resolve this inconsistency
in SMC averaging. Our experience shows that very good results can be obtained if only a small
portion of data before the time step at which the temperature change is applied (i.e. 400001), is
involved in SMC averaging. Here, pressure data of 30 time steps are enough. This means that
SMC averaging uses the data from the time step 399971 and after that. Convergence rates of
pressure in all of the bins are shown in Fig. 5. As can be seen, pressure increases from 2.9MPa
at temperature of 120K to 7.90MPa and 7.15MPa at temperatures of 180Kand 150K in bins 3



Analysis of pressure behavior in a temperature controlled... 887

Fig. 5. Convergence of pressure in different bins of the periodic flowwith the mean flow velocity
of 5m/s. A temperature gradient of 30K is enforced between bin 3 with 180K and bin 7 with 150K

at the time step 400001

and 7, respectively. It can also be seen that pressures in other bins are between these two values.
In the next Section, results of this case including pressure, density, and velocity gradients in
the periodic flow will be discussed in details. At this point, we conclude that pressure can be
effectively set to a desired value at anytime or anywhere in the solution domain by setting the
temperature of molecules. This finding can be very useful in applying flow boundary conditions
in the inlet or the outlet of the solution domain.

6. Creating a periodic flow with a temperature gradient

In this Section, we would like to analyze a flow field and gradients of its variables in a periodic
flow that is created by a temperature gradient. Consider the flowfield from the previous Section
that has been created by a temperature gradient enforced between bins 3 and 7. The mean
velocity was set to 5 m/s in the third bin, and temperature values of 180K and 150K were
set in bins 3 and 7, respectively. Having assumed periodic boundary conditions at the inlet and
outlet, bins 3 and 7 fall exactly in the middle of the flow field. It means that the distance from
bin 3 to bin 7 is equal to the distance from bin 7 to bin 3 in the periodic domain. From the time
step 400001, 800000 time steps of 10−15 s are taken to arrive at the equilibrium in the flow field.
The implemented temperature gradient of 30Khas resulted in pressure,mean flow velocity, and
molecular density variations shown in Fig. 6. In Fig. 6a, the mean flow velocity decreases from
5m/s in bin 3 to 3.3m/s in bin 7. Furthermore, pressure has reached its maximum value of
7.9MPa in bin 3 and its minimum value of 7.15MPa in bin 7 (Fig. 6b).

Fig. 6. Variations of flow variables along the periodic flow generated by temperature gradient of 30K,
(a) velocity, (b) pressure and (c) molecular density



888 H.R. Najafi, S.M.H. Karimian

As shown in Fig. 6c, the lowest value ofmolecular density occurs in bin 3where temperature
value is set to 180K. This is the highest temperature in the channel. Themolecules move with
a high virial velocity in this bin and, therefore, they push each other to the outside of the bin.
This is why the lowest molecular density occurs here. According to Eq. (3.1), while pressure
increases with temperature only through the first term and the effect of molecular density
appears in both terms, temperature has the dominant role in determining the pressure value
and its variation. Therefore, as seen in Fig. 6b, the highest pressure occurrs at the bin with the
highest temperature.

In the next step, we would like to investigate the effect of the temperature gradient on the
flow field at different temperature levels. We believe that the temperature level at which the
temperature gradient is applied, is important. Variations of velocity, pressure and molecular
density for a temperature gradient of 30K between bins 3 and 7, and at different fluid tempera-
ture values are shown in Fig. 7. The results are obtained for fixed temperatures of 150K, 180K,
210K, 240K, and 270K at bin 3. Themean flow velocity of 5m/s is still forced at bin 3.

Fig. 7. Variations of flow variables along the channel generated by a temperature gradient of 30K at
different fluid temperatures, (a) mean flow velocity, (b) pressure and (c) molecular density

In comparison to the level of fluid temperature, the results show that a 30K temperature
gradient becomes less effective at higher fluid temperature levels. In other words, the 30K
temperature gradient at the fluid temperature of 150Khas a higher effect on the flow characters
than it does at a fluid temperature of 270K, for instance. Therefore, it can be seen in Fig. 7
that as the temperature level increases, the mean flow velocity profile and molecular density
profile become smoother. In a similar fashion, the pressure profile becomes almost flat at higher
temperature values. Furthermore, in Figs. 7a and 7c, the profiles converge to a single avergage
value. This is not the case with the average pressure value which continuously increases with
temperature throughout the channel.

From Eq. (3.1), one can conclude that the mean flow velocity will not affect the pressure
values since it hasnodirect contribution in the equation.To study this, thevalues of temperature
in bins 3 and 7 are set to 150K and 300K, respectively enforcing a high temperature gradient
of 150K in the domain. The molecular flow is solved for this temperature gradient at different
mean flow velocities of 5, 10, 15, 20, 25, 30, 35, 40, and 45 meters per second in bin 3. The
profiles of velocity, pressure andmolecular density as well as temperature in these experiments
are shown in Figs. 8a-8d. In addition, the average temperature value of the molecules in the
whole domain has been calculated and plotted in Fig. 8e.

It can be seen in Fig. 8b that pressure decreases uniformly with an increase in the mean
flowvelocity. This is in contrastwith the above-mentioned statement. To further understand the
reason behind this behavior, let us analyze the results shown in Fig. 8. FromFigs. 8c and 8d, we
observe thatmolecular density is at its maximum in the lower temperature zone of the domain.



Analysis of pressure behavior in a temperature controlled... 889

Fig. 8. Variations of flow variables along the channel generated by a temperature gradient of 150K at
different mean flow velocities. Variations of (a) velocity, (b) pressure and (c) molecular density,

(d) temperature, (e) average temperature of molecules in the whole channel

This trend is more pronounced when the mean flow velocity increases. This means that the
number of molecules with lower temperature values increases in the solution domain with the
mean flow velocity. As a result, the bulk molecular temperature of the domain decreases with
the mean flow velocity (Fig. 8e). Since temperature has the dominant influence on pressure, it
can be justified that the reduction of pressure in the domain is a direct result of the reduction
in the bulkmolecular temperature.

The effect of the temperature gradient on the flow field is also studied here. Different tem-
perature gradient values are examined by setting temperature of bin 7 at 175K, 200K, 225K,
250K, 275K, 300K, 325K, and 350K. Amean flow velocity of 5m/s and temperature of 150K
are set at bin 3. Variations of velocity, pressure andmolecular density for different temperature
gradients between bins 3 and 7 are shown in Fig. 9.

At higher temperature values, the virial velocity of molecules increases and, therefore, mole-
cules repel each other. As a result, the molecular density decreases. In addition, the mean flow
velocity of molecules increases since at low molecular density values, fewer obstacles exist on
their ways. Based on this argument, the behavior of velocity and density profiles in Figs. 9a
and 9c can be easily justified.

As seen in Fig. 9b, pressure in the whole domain rises with the temperature gradient. An
increase in the temperature gradient increases the average temperature of the flow field globally
and, therefore, pressure rises in the whole channel. The behavior of pressure profiles can be
explained locally as well. In bin 3 for instance, pressure rises due to an increase in molecular
density although temperature is constant. In bin 7, however, pressure rises due to an increase
in temperature althoughmolecular density reduces. Once again, note that temperature has the



890 H.R. Najafi, S.M.H. Karimian

Fig. 9. Variations of flow variables along the channel generated by different temperature gradients at a
mean flow velocity of 5m/s in bin 3. Variations of (a) velocity, (b) pressure and (c) molecular density

primary role in variation of the pressure value and its behavior. Therefore, pressure increases
from bin 3 to bin 7.

The last test case has been repeated with a fixed temperature of 180K (instead of 150K) to
investigate the flow behavior at a different temperature level. The results are not reported here
to limit the length of the present paper; but similar behavior has been observed.At the end, it is
worth to mention that while more accurate methods for pressure calculation in inhomogeneous
fluids can be found in the literature (Todd et al., 1995), we believe that employment of such
methods will not lead to different conclusions than those made herein.

7. Conclusion

Details of pressure calculation in amolecular system and its sensitivity to approximations have
been studied herein. It has been shown that for molecules near boundaries, inclusion of all
molecules within their cut-off region including those outside of the domain boundary is crucial
for correct calculation of pressure.

In this paper, a method has been introduced to create a flow by controlling macroscopic
properties in twodifferent regions.Temperaturevalueshavebeen imposedat two selected regions
along the flow, and the mean flow velocity has been imposed at one of these regions. Having
applied these conditions to theperiodicflowandconductingparametric studyon it, the following
results have been obtained.

• Temperature and its gradient have themost dominant role on the variation of the pressure
value and its gradient in the periodic flow.

• In the periodic flow, the pressure gradient established as a result of a constant tempe-
rature gradient is not always the same at different temperature levels. In fact, both the
temperature gradient and the temperature atwhich this gradient is applied, determine the
pressure gradient in the flow.

• Since the mean velocity does not appear in the pressure formula, it is expected that
variations of the mean velocity would not change pressure in the periodic flow. However,
the results show that for a constant temperature gradient, pressure changes inversely with
the mean flow velocity. Based on our analysis, this is because the mean flow velocity
directly changes the average temperature of the flow within the channel.



Analysis of pressure behavior in a temperature controlled... 891

References

1. Allen M.P., Tildesley D.J., 1987,Computer Simulation of Liquids Oxford, New York

2. Branam R.D., Micci M.M., 2009, Comparison of wall models for the molecular dynamics simu-
lation of microflows,Nanoscale and Microscale Thermophysical Engineering, 13, 1-12

3. Darbandi M., Abbassi H.R., Khaledi Alidusti R., Sabouri M., 2011,Molecular dynamics
simulation of nano channel as nanopumps, ICNMM, Edmonton, Alberta, Canada

4. Fan X.J., Phan-Thien N., Teng Yong N., Diao X., 2002,Molecular dynamics simulation of
a liquid in a complex nano channel flow,Physics of Fluids, 14, 3, 1146-1153

5. HanM., 2008,Thermally-drivennanoscalepumpbymoleculardynamics simulation,JournalTitle:
Journal of Mechanical Science and Technology, 22, 157-165

6. Hanasaki I., Nakatani A., 2006, Fluidized piston model for molecular dynamics simulations of
hydrodynamics flow,Modelling and Simulation in Materials Science and Engineering, 14, S9-S20

7. HuangC., Choi P.Y.K., NandakumarK., Kostiuk L.W., 2006,Molecular dynamics simula-
tion of a pressure-driven liquid transport process in a cylindrical nanopore using two self-adjusting
plates, Journal of Chemical Physics, 124, 234701

8. HuangC., NandakumarK., Kwok D.Y., 2004,Non-equilibrium injection flow in a nanometer
capillary channel, ICMENS’04, 374-378

9. Kamali R., Kharazmi A., 2011, Molecular dynamics simulation of surface roughness effects on
nanoscale flows, International Journal of Thermal Sciences, 50, 3, 226-232

10. Karniadakis G.E., Beskok A., Aluru N., 2002, Micro Flows and Nanoflows, Springer, New
York, 641-648

11. Karimian S.M.H., Izadi S., 2013, Bin size determination for the measurement of mean flow
velocity inmolecular dynamics simulation, International Journal ForNumericalMethods In Fluids,
71, 7, 930-938

12. Karimian S.M.H., Izadi S., Barati Farimani A., 2011, A study on themeasurement of mean
velocity and its convergence in molecular dynamics simulations, International Journal for Nume-
rical Methods in Fluids, 67, 12, 2130-2140

13. Karimian S.M.H., Namvar S., 2012, Implementation of SMC averagingmethod in a channeled
molecular flow of liquids and gases, Journal of Physics: Conference Series, 362, 1, 2029

14. Kim B.H., Beskok A., Cagin T., 2010, Viscous heating in nanoscale shear driven liquid flows,
Microfluid Nanofluid, 9, 31-40

15. Koplik J., Banavar J.,Willemsen J., 1988,Molecular dynamics of poiseuille flow andmoving
contact lines,Physical Review Letters, 60, 1282-1285

16. Koplik J., Banavar J. R., Willemsen J.F., 1989, Molecular dynamics of fluid flow at solid
surfaces,Physics of Fluids A, 1, 781-794

17. Leach A.R., 1999,Molecular Modeling: Principles and Applications, Longman

18. Liu C., Li Z., 2010, Molecular dynamics simulation of composite nanochannels as nanopumps
driven by symmetric temperature gradients,Physical Review Letters, 105, 174501

19. Mi X.B., Chwang A. T., 2003, Molecular dynamics simulations of nanochannel flows at low
Reynolds numbers,Molecules, 8, 193-206

20. Nagayama G., Cheng P., 2004, Effects of interface wettability onmicroscale flow bymolecular
dynamics simulation, International Journal of Heat and Mass Transfer, 47, 501-513

21. Namvar S., Karimian S.M.H., 2012, Detailed investigation on the effect of wall spring stiffness
on velocity profile in molecular dynamics simulation, Journal of Physics: Conference Series, 362,
1, 2039



892 H.R. Najafi, S.M.H. Karimian

22. Plimpton S.J., 1995, Fast parallel algorithms for short-range molecular dynamics, Journal of
Computational Physics, 117, 1

23. Priezjev N.V., 2007, Effect of surface roughness on rate-dependent slip in simple fluids, Journal
of Chemical Physics, 127, 144708

24. Rapaport D.C., 2004,The Art of Molecular Dynamics Simulation, Cambridge University Press

25. Sadus R.J., 2002, Molecular Simulation of Fluids: Theory, Algorithms and Object-Orientation,
Elsevier

26. Sofos F., Karakasidas T.E., LiakopoulosA., 2009,Non-equilibriummolecular dynamics in-
vestigation of parameters affecting planar nanochannel flows,Contemporary Engineering Sciences,
2, 6, 283-298

27. StillingerF.H.,RahmanA., 1974, Improved simulation of liquidwater bymolecular dynamics,
Journal of Chemical Physics, 60, 1545-1557

28. Sun M., Ebner C., 1992, Molecular-Dynamics simulation of compressible fluid flow in two-
-dimensional channels,Physical Review A, 46, 4813

29. Todd B.D., Evans D.J., Davis P.J., 1995, Pressure tensor for inhomogeneous fluids, Physical
Review E, 52, 1627

30. Travis K.P., Evans D.J., 1997, Molecular spin in a fluid undergoing Poiseuille flow, Physical
Review E, 55, 1566-1572

31. Travis K.P., Gubbins K.E., 2000, Poiseuille flow of Lennard-Jones fluids in narrow slit pores,
Journal of Chemical Physics, 112, 1984-1994

32. Tysanner M.W., Garcia A.L., 2004, Measurement bias of fluid velocity in molecular simula-
tions, Journal of Computational Physics, 196, 173-183

33. Tysanner M.W., Garcia A.L., 2005, Non-equilibrium behavior of equilibrium reservoirs in
molecular simulations, International Journal of Numerical Methods in Fluids, 2050, 1-12

34. Verlet L., 1967, Computer “experiments” on classical fluids. I. Thermodynamical properties of
Lennard-Jonesmolecules,Physical Review, 159, 98

35. Zhang Z.Q., Zhang H.W., Ye H.F., 2009, Pressure-driven flow in parallel-plate nanochannels,
Applied Physics Letters, 95, 154101

36. Ziarani A.S., Mohamad A.A., 2005, A Molecular dynamics study of perturbed Poiseuille flow
in a nanochannel,Microfluid Nanofluid, 2, 12-20

Manuscript received February 12, 2015; accepted for print December 10, 2015