Thermodynamic characterization of two layers of CO2 on a graphite surface

(1)

ContentslistsavailableatScienceDirect

Chemical

Physics

Letters

j o ur na l ho me p ag e :w w w . e l s e v i e r . c o m / l o c a t e / c p l e t t

Thermodynamic

characterization

of

two

layers

of

CO

2

on

a

graphite

surface

T.T.

Trinh

a

_,

_D.

_Bedeaux

a

_,

_J.-M.

_Simon

b

_,

_S.

_Kjelstrup

a,c,∗

a_Department_of_Chemistry,_Norwegian_University_of_Science_and_Technology,_Trondheim,_Norway

b_Laboratoire_{Interdisciplinaire}_Carnot_de_Bourgogne,_UMR-6303_{CNRS-Université}_de_Bourgogne,_Dijon,_France

c_Department_of_Process_and_Energy_Laboratory,_Delft_University_of_Technology,_Delft,_The_Netherlands

a

r

t

i

c

l

e

i

n

f

o

Articlehistory:

Received3July2014

Inﬁnalform8August2014

Availableonline14August2014

a

b

s

t

r

a

c

t

Wefindbyexaminationofdensityprofilesthatcarbondioxideadsorbsongraphiteintwodistinctlayers. Wereporttheactivitycoefficient,entropyandenthalpyforCO2ineachlayerusingaconvenient

com-putationalmethod,theSmallSystemMethod,therebyextendingthismethodtosurfaces.Thisopensup thepossibilitytostudythermodynamicpropertiesforawiderangeofsurfacephenomena.

1. Introduction

Phasetransformations[1],formationsofnanostructures[1,2], metastablephases[3]oragglomerations[4],aswellaschemical reactions[5]areallphenomenawhichbecomespecialatsurfaces. Asurfacecanalsoposeamajorbarriertotransport[6,7]. Thermody-namicdatawillimprovetheunderstandinganddynamicmodeling ofequilibriumandnonequilibriumstates.Butthermodynamicdata forsurfacesarenoteasilyavailablebyexperiments.Withthe intro-ductionoftheSmallSystemMethod[8,9]computationalresults havebecomefeasible.Thismethodhassofarnotbeenusedto studysurfaces.Inthisworkweextendthemethodtoatechnically importantprocess,thephysisorptionofCO2onagraphitesurface.

Graphiticmembranesarepromisingcheapcandidatemembranes forCO2separationpurposes[10].Theprocesscanbewritten:

CO2(g)+graphiteCO2(s) (1)

where (g) means gas phase and (s) means adsorbed gas. The purposeoftheworkistodemonstratehowthermodynamic infor-mationcanbegainedofanadsorbedstateusingthenewsimulation method.Weshallseethatadsorptiontakesplaceintwodistinct thermodynamiclayers,andthattheﬁllingintothelayerscanbe regardedasatrade-offbetweentheentropiesandenthalpiesfound foreachlayer.

∗ Correspondingauthorat:DepartmentofChemistry,NorwegianUniversityof

ScienceandTechnology,Trondheim,Norway.

E-mailaddresses:trinhthanhthuat@gmail.com(T.T.Trinh),

signe.kjelstrup@ntnu.no(S.Kjelstrup).

2. Methods

The Small System Method [8,9] makes use of the relation betweentheinverseofthermodynamiccorrectionfactorandthe ﬂuctuationinparticlenumberN.Forasurface,therelationis

1 =

N2

₋

_N

2

N

T,,A (2) Fluctuationsaresampledinsmallopensystems(disks) with areaAinsideareservoir.Thereservoiriscreatedasalarge rect-angularboxwithperiodicboundaryconditions.Thetemperature, T,andthechemicalpotential,,inthereservoir arecontrolled. Thethicknessofthesamplingsystemissettothethicknessofthe surface(seebelow).Theareaofthesamplingsystemisvaried, vary-ingtheradiusofthedisk,L.Thesmallestradiususedissosmall thatitallowsonlyonemoleculeinsidethedisk;thelargestradius contains20–30molecules.Theinversethermodynamicfactorisa linearfunctionofthe1/L[11–13]inaparticularrange,tobefound foreachcase:

1 = 1 ∞,s

1+B L

(3) HereBisasmallsystemspeciﬁcconstantandsuperscript∞means thevalueinthethermodynamiclimit.Byextrapolatingthelinear regimetothethermodynamiclimit,weobtainthewantedquantity ∞,s.Lateronweuses₌∞,s_for_simplicity.

Thechemicalpotentialofagasadsorbedtoasurfaceinalayer iis: s i= 0,s i +RTlna s i= 0,s i +RTln s i Cs i C0,s_i (4) http://dx.doi.org/10.1016/j.cplett.2014.08.026

(2)

Here 0,s_i is thestandard chemical potential, as _is _the

(dimen-sionless)activityoftheadsorbedphase,Cs _is_the_surface_excess

concentration,ands_is_the_activity_{coefﬁcient.}_The_standard_state

is thehypothetical idealstate havings₌₁ _at_layer _saturation,

C0,s _(given_in_particles_per_(nm2_))._The_entropy_follows_from_the

standardrelationSs₌₋_(ds_/dT)_and_the_enthalpy_of_the_layer_from

Hs i=

s i+TS

s

i. Weshallﬁndtheseproperties fortheadsorption

(1).Totalsurfaceexcessconcentrations(adsorptions)aredeﬁned accordingtoGibbs[14,15]

Cs₌

∞

0

(C(z)−Cg_(˛)(z₋_˛))_dz ₍₅₎

whereC(z)andCg_are_the_{concentrations}_of_adsorbed_molecules

andofmoleculesinthegasphase,respectively.Theintegrationis carriedoutfromtheequimolarsurfaceofgraphite(atz=0)tothe bulkgasphase(z=∞).TheHeavisidefunction,,isbydefinition unitywhentheargumentispositive,andzerowhentheargument isnegative.Fordefinitionoflayerconcentrations,seeSection3.The thermodynamicfactorisalsodefinedby:

s=1+

∂

lns

∂

lnCs

T,A (6) WecandeterminetheactivitycoefﬁcientbyintegratingEq.(6) fromzeroadsorption,oncethethermodynamicfactorisknown fromtheSmallSystemMethod.Thisgives

s 1 dlns₌

Cs 0 (s₋_1)d_ln_Cs ₍₇₎

WhenCs_→_0,_there_is_no_interaction_between_particles,_meaning

thats₌_1._We_also_have

s−s=

s s ds=RT

Cs Cs sdlnCs (8)

meaningthatwecanﬁndthechemicalpotentialatany adsorp-tionrelativetoareferencestate(indicatedby)fromthisequation. Weshallusethisrelationtodeterminethestandardstate chem-icalpotential0,s_i byplottingthelefthandsideofEq.(8)versus ln(s

iC s i/C

0,s

i ),extrapolatingtothestatewherethisexpressionis

zero.

3. Simulationdetails

ThesystemconsistedofsheetsofgraphiteandCO2molecules.

Thegraphitewasmadefrom5sheetsof graphenewithoutany defects.We orientedthesheets in theboxsuchthat thesheet surfaceswereperpendiculartothez-axis.Thedistancefromthe graphitesurfacewasmeasuredalongthisaxis,takingthe equimo-larsurfaceofgraphiteaszero.Thegraphitelayerswereﬁxedin space,stillyieldinggoodresultsforadsorptionanddiffusionofgas onthesurface[16,17].Arigidbodymodel,TraPPE,wasusedfor CO2 [18].TheintermolecularpotentialbetweenCO2–CO2 wasa

shiftedandtruncated12-6Lennard–Jones(LJ)potential[19]with long-rangeCoulombinteractions,whichweredealtwithusingthe Ewaldsummationtechnique[19].TheinteractionbetweenCO2and

theC-atomsofgraphitewassimilarlydescribedbyaLJpotential. Detailsofparametersforthesimulationweregivenearlier[16]. Theyhavebeenconﬁrmedtoyieldanaccurateadsorptionenergy forCO2ongraphitesurface[16].

Classicalmoleculardynamics(MD)simulationswereperformed usingtheLAMMPSpackage[20].Asnapshotofthesystemisgiven in Figure1,showing graphiteat 500K withCO2 ina relatively

densegasandinanadsorbedstate.Periodicboundaryconditions wereusedforthereservoirinalldirections.Thesimulationwas donewithtime steps of 0.001ps. The initialconﬁgurationwas

Fig.1.SnapshotshowingCO2inthegasphaseataﬂuid-likedensityandadsorbed

onagraphitesurfaceat500K.ThenumberofCO2particlesinthesystemwasNCO2=

2800.Thegreenandredcolorsrepresentcarbonandoxygenatoms,respectively.

(Forinterpretationofthereferencestocolorinthisﬁgurelegend,thereaderis

referredtothewebversionofthisarticle.)

constructedby randomlydistributing CO2 molecules above the

graphitesurface.Thesystemwasstabilizedduring2000psbyruns withconstantNtotVtotTusingNosé–Hooverthermostats[21].

InordertochecktherangeofvalidityofEq.(3),thereservoir sizewasvaried.Thereservoirhastobelargeenoughcomparedto thesamplingsystem.Forahardspheresystem,wefoundthata suitablereservoirhadalength20timesthehardsphereradiusL [12,13].Simulationboxsizesof(a=42 ˚A,b=54 ˚A,c=84 ˚A),(2a,2b, c),and(3a,3b,c)weretested.Thelargestsystemsgaveidentical results,sobox(2a,2b,c)wasused.

Weperformed2000psrunsinareservoiratthermal equilib-rium,varyingLfrom1.3to20 ˚A,samplingrandomly30disksacross thegraphitesurface.Thetotalnumberofframeswas10000.Hence, thevalueofs_for_each_L_was_obtained_from_the_statistics_of₃₀₀₀₀₀

samples[11,13].Thesamplingwasdonefortemperaturesranging from300Kto550K,andforanumberofCO2particlesinthebox,

varyingfrom200to2800(correspondingtoapressurerangeof 1–60barat300K).

4. Resultsanddiscussion

4.1. Adsorptions

Wedescribeﬁrsthowweﬁndtheadsorptions Cs

1 andC2s for

eachlayer.ThetotaladsorptionisthenCs₌_Cs

1+C2s.The

adsorp-tion(excessconcentration)ofCO2wasfoundusingtheresultsin

Figure2withthedeﬁnition(5).Theﬁgureshowsthedensityof CO2 moleculesasafunctionofthedistancetothegraphite

sur-face,andisaquantiﬁcationofsnapshotsliketheoneshowedin Figure1.WeobserveinbothﬁguresthatCO2canformtwolayers

onthegraphiticsurface.InFigure2,weseealayeraroundthefirst peakwhichextendsfromthegraphitesurfacetoposition␣,while asecondpeakextendsfrom␣to␤.WerefertotheintegralinEq. (5)duetothefirstpeakasthebottomlayeradsorption,whilethe correspondingintegralrelatedtothesecondpeak,iscalledthetop layeradsorptiononthegraphite.Thepositionofthefirstpeakat 5 ˚Acoversapproximatelyonemolecularlayer,asthelengthofthe moleculeis5.4 ˚A.Thesnapshotshowsthatmostofthecarbon diox-idemoleculesarelyingparalleltothesurface.Thisgivesaheight morelikethediameterofanoxygenion,3.6 ˚A.Athicknessof5 ˚A canbeinterpretedtoincludesomemoleculeswhicharestanding slightlytiltedwithrespecttothesurface,andsuchmoleculescan alsobefoundbyvisualinspectionofFigure1.Itisinterestingthat thepositionofplane␣isalwaysnear5 ˚A,independentofthe tem-perature.Thiscanreflectthatmoleculesareeitherlyingorstanding

(3)

Fig.2.ThedensityofCO2moleculesasafunctionofdistancetothesurfaceina

reservoirwithNCO2=2800.Thetemperaturesare500Kand350K.Wedistinguish

betweenthreezones,from0to␣:ﬁrstadsorbedlayer,␣to␤:topadsorbedlayer,

above␤:gasphase.Thebottomlayerextendstoaround5 ˚A.Thethicknessofthe

toplayerislargerat350Kthanthatat500K.

inthefirstlayer.Theattractiveforcesofthegraphite,notbeingable toreachabovethelayer,seemtobecentralforthislayer.Butthe factthatthetoplayerstartstobefilledbeforethebottomlayeris full,motivatesadivisionofthewholesurfaceintotwolayers.The toplayer,whichextendsfrom␣to␤,appearsalsoratherdifferent fromthefirst.WeseefromFigure2thatthepositionoftheplane␤ varieswithtemperature,unlikethepositionofplane␣.Forthe sys-temwith2800CO2molecules,␤issmallerat500Kthanat350K.

Thetoplayeristhusmorediffuse.Fewattractiveforcesareableto keepthemoleculeswithinthissurfacelayer,whenthemolecular kineticenergybecomeslarger.TheadsorptionofCO2ongraphiteis

alwayshigheratlowtemperatureinbothlayers,meaningthatthis kindoftrade-offbetweenkineticandpotentialenergiesappliesto bothlayers.

Weﬁndthatthebottomlayerhasalwaysmore thanhalfof thetotal adsorption, independentof temperature (not shown). Thefractionislarger,thelowerthetemperatureis,butdoesnot exceed0.9.Thefractiondecreaseswithincreasingtotaladsorption, becausethebottomlayerbecomessaturated,whilethetoplayer keepsgrowing.Fullcoverage(maximumadsorption)ofthe bot-tomlayerisobtainedat300KwithC₁0,s=12.5(molecules/nm2_).

Wethereforechoseasstandardstate,theidealstatewithC₁0,s= 12.5(molecules/nm2_), _{corresponding} _to_0.31 _CO

2 moleculeper

surfacecarbonatom.

4.2. Thethermodynamiccorrectionfactor

Resultsfortheinversethermodynamicfactor(s₎−1_were

plot-tedversus 1/L for the top and bottom layers separateand for theircombination (notshown).Thecurvesapproached unityas expected,when1/Lwasincreased.Thisisthesmallsystemlimit, whichhasatmostoneparticleinthesamplingsystem.Inorderto ﬁndthethermodynamiclimitvalue,weappliedlinearregression topointsintheinterval0.1<1/L<0.4andextrapolatedthelineto largeL.Theregionusedforextrapolationcoincidedwiththeregion foundearlier[11,13].

Thethermodynamiclimitvalueoftheinversethermodynamic correctionfactor(s₎−1_was_plotted_as_a_function_of_the_adsorption

ineachsinglelayer,andinthecombinedlayers,forall temper-atures.AtypicalexampleisshowninFigure3forT=500K.The bottomlayerhadalwaysthesmallest(s₎−1_,_while_the_total_layer

hadthebiggestvalueof(s₎−1_._We_found_that₍s₎−1_decreased

more or less linearly with the adsorption. The results for the

Fig.3.Thermodynamiclimitvaluesfor(s₎−1_in_the_top,_bottom_and_total_CO2_layer

ongraphiteat500K.Astraightlineisﬁttedtotheresultsandforcedthrough1on

they-axis.

bottomlayerandthetwolayerscombined(thetotal layer) fol-lowedthestraightlinenicely,whilethescatterofresultsaround thelinedrawnforthetoplayerwaslarger.

Allcurvesmustextrapolatetoathermodynamiccorrection fac-torequal1,whentheadsorptiongoestozero.Thiswasfoundthe relationbetween(s₎−1_of_the_total_layer_and_the_total_adsorption

wasquitelinear.Theobservedlineardependenciesmeanthatthe inversethermodynamicfactorofthetotallayer,canbefoundby addingcontributionsfromthebottomandtoplayers.

4.3. Theactivitycoefficientsforcarbondioxideinthesurfaces Asetoftypicalactivitycoefficients,obtainedat500KfromEq. (7),areshowninFigure4.Thecoefficientsapproach1asexpected atlowdensity,somewithmorenoisethanothers.Theincreasein thecoefficientforthebottomlayerislargerthanforthetoplayer andthereforelargerthanforthetotallayer.Thisreflectsthatthe moleculesarefurtherapartinthetoplayerthaninthefirst,making repulsiveforceslessrelevantinthetopthaninthebottomlayer. Themoleculesdonotadsorb,unless thesurfacebindingenergy canovercometherepulsion,however.Thetotallayeractivity coef-ficientwasexpressedasafunctionoftotalsurfaceadsorptionat differenttemperatures,witharegressioncoefficient0.99orbetter: s

i =1+aCis+b(Cis) 2

(9)

Fig.4.Activitycoefﬁcientofbottom,topandtotallayerasafunctionofthelayer

(4)

Table1

ParametersthatdescribethetotallayeractivitycoefﬁcientintheempiricalEq.(9).

Theregressioncoefﬁcientwas0.99orbetter.

T(K) b a 300 0.0060 0.0005 350 0.0059 0.0054 400 0.0055 0.0066 450 0.0059 0.0233 500 0.0057 0.0320 550 0.0067 0.0216

Fig.5. ChemicalpotentialofCO2onagraphitesurfaceasafunctionofln(sCs/C0,s)

atdifferenttemperatures.ThestraightlineisﬁttedtoEq.(4).

ParametersaandbaregiveninTable1.Theseresultsmaybe usefulforthermodynamicmodelingofthetotalCO2adsorbedona

graphitesurface(Fig.5).

4.4. Thechemicalpotentialofthesurfaceadsorbedgas

Thechemical potentialof theadsorbed staterelative tothe chemicalpotentialatthelowestadsorptionwasfoundfromEq. (8). The chemical potential difference for the total layer was plotted vs the right hand side for the different temperatures used.UsingtherelevantactivitycoefﬁcientandC₁0,s=C₂0,s=C0,s₌

12.5molecules/(nm2_), _we _next _used _Eq. ₍₄₎ _to _plot _the _same

data.Theslope ofthelinear ﬁt gaveRTas expectedwithgood accuracy.Thestandardstatevalue,0,s_,_was_found_for_each

tem-peraturefromthiscurvesettinglns i(C

s i/C

0,s

i )=0. Forthetotal

layer,weobtained0,s_(T0₎₌_−7.6_kJ/mol_at₂₉₈_K._The_temperature

variationofthechemicalpotentialwasusedtoﬁndthestandard entropy S0,s₌_10.6_J(mol_K). _It _was _constant _in _the _temperature

intervalused.Thestandardenthalpyat298Kwasthencalculated toH0,s₌_−4.4_kJ/mol._Data_for_the_total_and_separate_layers_are_given

inTable2.

Inordertocheckthattheresultsfortheseparatelayerswere internallyconsistent,weexpressedthechemicalpotentialsofthe topandbottomlayersofCO2 byEq.(2).Atequilibrium,wehave

s

1=s2,wherethesubscriptdenoteslayernumber.Byintroducing

expression(2)intothisequilibriumcondition,weobtained

Cs 1 Cs 2

=

s 2 s 1

e(0,s2 − 0,s 1 )/RT ₍₁₀₎

By plotting the left hand side versus the ratio of activity coefﬁcients,ﬁttingtheplottoastraightline,wecalculated0,s₂ − 0,s₁ fromtheslopeasafunctionoftemperature.Thedifferencein thestandardchemicalpotentialswasconsistentwiththedifference obtainedfromthedatainTable2.

Thethermodynamicdatathatwehavedeterminedforthetotal layerofadsorbedgas,arereasonable.Thestandardchemical poten-tialandtheenthalpyarebothnegative,indicatingthatadsorption is favorable in the standard state,in spite of a largereduction in the entropy from the gas phase tothe surface, a reduction whichislargerthanthatfromthegastotheliquidstate.Wehave furthermoreseenthatwecandistinguishbetweentwoseparate thermodynamically deﬁned layersof carbondioxidewithin the wholelayeronthegraphitesurface.Eachlayerhasitsown ther-modynamicproperties.WeseefromTable2,thatthebottomlayer entropyisrelativelysmallcompared tothegasentropyofCO2.

Infactitcancomparetothecarbonentropyinagraphitelattice. ThenegativeentropydifferencerepresentedbytheorderingofCO2

whenitisadsorbedmustbeovercomebyarelativelylargeenthalpy foradsorptionintothelayer.Thetoplayerentropyislargerthan theentropyofthebottomlayer.Itisnotgas-like,butfourtimes higherthanthatofthebottomlayervalue.Thislayercantherefore bestabilizedwithasmallerenthalpy.Thetotallayerpropertiescan besaidtomaskthepropertiesofthesinglelayers,asthetotallayer hasstandardentropyclosertothetoplayer,andanenthalpymore likethebottomlayer.Informationoftheseparatelayersobtained fromthesimulationwillthereforeaddinsighttothesystem.Such insightcannotbeeasilyobtainedbyexperiments,asitisdifﬁcult todistinguishbetweenthelayersinanadsorptionexperiment.It ishoweverrathersimpletocompareanexperimentallyobtained adsorptionisothermtotheadsorptionisothermderivedherefor thetotaladsorption.

5. Conclusion

WehaveshownthatarecentlydevelopedSmallSystemMethod canbeusedtocalculateaconsistentsetofthermodynamicdatafor surfacesfromoneMolecularDynamicssimulation.Forcarbon diox-ideadsorbedtographite,wedeterminedthechemicalpotential, theactivitycoefﬁcient,theentropyandenthalpydirectly.Closer inspectionofdensityproﬁlesrevealedthatwecanspeakoftwo, notonelayerofadsorbedgas.Thethermodynamicanalysesreveal adistinctionbetweenthelayers.Thegashaslargerentropyinthe topthaninthebottomlayer,ismoremobileinthislayer.Thetwo layersareinequilibriumwithoneanother,meaningthata lower-ingoftheentropycantakeplace,iftheenthalpycancompensate forthechangeinentropywhenatopparticlemovestothebottom layer.Suchinformationisinvaluableinthemodelingand expla-nationofsurfaceprocesses.Extensionstoothersurfacesshouldbe straightforward.

Table2

ThermodynamicdataforstandardstateadsorptionofcarbondioxideatagraphitesurfaceatT=298K.Theuncertaintyinthedeterminationsisontheaverage3%.

Correspondingvaluesforgraphiteandcarbondioxidegasarealsogiven.

T=298K Totallayer(thiswork) Bottomlayer(thiswork) Toplayer(thiswork) Graphitea _CO2_gasa

0,s i (kJ/mol) −7.6 −9.8 −5.1 – – S0,s i (J/molK) 10.6 4.1 13.7 5.74 213.78 H0,s i (kJ/mol) −4.4 −8.6 −1.0 1.05 9.36

(5)

Acknowledgments

TheauthorsacknowledgeTheResearchCouncilofNorwayRCN projectno209337andTheFacultyofNaturalScienceand Technol-ogy,NorwegianUniversityofScienceandTechnology(NTNU)for ﬁnancialsupport.ThecalculationpowerisgrantedbyThe Norwe-gianMetacenterforComputationalScience(NOTUR).ETHZürichis thankedforguestprofessorshipstoSKandDB.

References

[1]L.P.H.Jeurgens,Z.M.Wang,E.J.Mittemeijer,Int.J.Mater.Res.100(2009)1281.

[2]F.Costache,S.Kouteva-Arguirova,J.Reif,Appl.Phys.A79(2004)1429.

[3]R.BradleyShumbera,H.H.Kan,J.F.Weaver,Surf.Sci.601(2007)235.

[4]D.J.Srolovitz,M.G.Goldiner,JOM:J.Min.Met.Mater.Soc.47(1995)31.

[5]M.P.Cox,G.Ertl,R.Imbihl,J.Rustig,Surf.Sci.134(1983)L517.

[6]I.Inzoli,J.M.Simon,S.Kjelstrup,Langmuir25(2009)1518.

[7]I.Inzoli,J.M.Simon,D.Bedeaux,S.Kjelstrup,J.Phys.Chem.B112(2008)14937.

[8]X.Liu,S.K.Schnell,J.M.Simon,D.Bedeaux,S.Kjelstrup,A.Bardow,T.J.Vlugt,J. Phys.Chem.B115(2011)12921.

[9]P.Krüger,S.K.Schnell,D.Bedeaux,S.Kjelstrup,T.J.Vlugt,J.-M.Simon,J.Phys. Chem.Lett.4(2012)235.

[10]X.He,J.ArvidLie,E.Sheridan,M.-B.Hägg,EnergyProcedia1(2009)261.

[11]S.K.Schnell,X.Liu,J.-M.Simon,A.Bardow,D.Bedeaux,T.J.H.Vlugt,S.Kjelstrup, J.Phys.Chem.B115(2011)10911.

[12]S.K.Schnell,T.J.H.Vlugt,J.-M.Simon,D.Bedeaux,S.Kjelstrup,Chem.Phys.Lett. 504(2011)199.

[13]S.K.Schnell,T.J.H.Vlugt,J.-M.Simon,D.Bedeaux,S.Kjelstrup,Mol.Phys.110 (2011)1069.

[14]J.W.Gibbs,TheScientiﬁcPapersofJWGibbs,vol.1,Dover,NewYork,1961.

[15]J.W. Gibbs,OntheEquilibriumofHeterogeneousSubstances,Connecticut AcademyofArtsandSciences,1877.

[16]T.T.Trinh,T.J.Vlugt,M.B.Hagg,D.Bedeaux,S.Kjelstrup,Front.Chem.1(2013) 38.

[17]T.S.vanErp,T.T.Trinh,S.Kjelstrup,K.Glavatskiy,Front.Phys.1(2013)30.

[18]J.J.Potoff,J.I.Siepmann,AIChEJ.47(2001)1676.

[19]M.P.Allen,D.J.Tildesley,ComputerSimulationofLiquids,OxfordUniversity Press,1989.

[20]S.Plimpton,P.Crozier,A.Thompson,LAMMPS–Large-scaleAtomic/Molecular MassivelyParallelSimulator,SandiaNationalLaboratories,2007.

[21]G.J.Martyna,M.L.Klein,M.Tuckerman,J.Chem.Phys.97(1992)2635.

[22]J.Cox,D.D.Wagman,V.A.Medvedev,CODATAKeyValuesforThermodynamics, Chem/Mats-Sci/E,1989.