Properties of the hard-sphere fluid at a planar wall using virial series and molecular-dynamics simulation

Paganini, Iván E.; Davidchack, Ruslan L.; Laird, Brian B.; Urrutia, Ignacio

doi:10.1063/1.5025332

我们研究了与平面硬壁接触的硬球流体。通过将非均匀维里级数与模拟结果相结合，我们获得了计算表面热力学性质（如表面吸附Γ和表面自由能（或表面张力））的新的精度基准，γ。我们简要介绍了选择分割面位置的问题，并通过建议使用Γ和Γ的替代函数来避免此问题γ与采用的参照系无关。最后，我们给出了系统表面热力学性质对堆积分数依赖性的解析表达式，确保了参数化函数在任何参考系下的高精度。Γ和Γ的拟议参数表达式γ，在统计误差范围内拟合准确的仿真结果。

I.简介

硬球（HS）模型是粒子间相互作用的最简单模型之一。在过去的一个世纪里，它在凝聚态物质物理概念的发展中发挥了关键作用。基于HS模型开展了大量研究：气体、液体和固体建模，以及冷凝、堵塞、玻璃化转变、结晶、固体缺陷和混合等研究。此外，HS模型有助于现代流体理论的发展。特别是，经常使用HS作为参考系统来研究非均匀流体，例如，使用经典密度泛函理论（DFT）的形式。¹

几十年来，人们一直致力于获得HS流体状态方程（EOS）的精确解析表达式，该表达式描述了体积压力对密度的依赖性。^2–11在积分方程的背景下，Percus-Yevick（PY）方程的解构成了这一方向的关键贡献。^12,13EOS的不同方法包括评估维里系数，以获得低密度下的分析相关性，并使用分子动力学（MD）模拟进行模拟^9,10,14–17和蒙特卡罗（MC）技术。¹⁸

除了对确定体积特性感兴趣外，随着时间的推移，密闭HS流体的特定特性也越来越受到关注。在这种情况下，焦点转移到壁流表面自由能EOS，γ和吸附Γ。维里级数的精神病学研究γ从贝勒曼开始，由索科·奥斯基和斯特基继续。^19–22最近，Yang、Schultz、Errington和Kofke（YSEK）对高阶维里系数进行了数值评估。¹⁸的行为γDavidchack、Laird和Roth（DLR）使用MD模拟研究了整个密度范围内直至冻结的情况。^23–27通过不同的技术获得的结果，如维里级数、MD、MC、DFT和标度粒子理论（SPT），允许发展对γ密度。^25,28–31然而，最近关于Γ和γ的行为的结果存在一些差异γ低密度（使用维里级数研究）和高密度（基于MD和DFT的不同研究之间存在差异）。

在这项工作中，我们对HS流体在平面硬壁上的表面自由能和吸附进行了新的精确计算。我们使用了一种双重方法来实现这一点，该方法结合了簇积分和维里系数的精确评估以及新的高分辨率MD模拟结果。这些技术是互补的，允许我们提出准确描述γ密度上的Γ。结果包括高达七阶的团簇积分和维里系数，以及HS流体整个密度范围内的模拟数据。以秒为单位。二，我们描述了非均匀平面系统的基本热力学方程，以及几种性质的幂级数（virial）表示。章节三介绍了用于计算簇积分和维里系数的方法，见第节。四、，我们描述了通过模拟获得表面热力学性质数据的方法γ和Γ。我们在第。V（V）其中，我们重点研究与分割表面位置无关的表面热力学性质，这使我们能够获得表面自由能和吸附的精确参数表达式。当分割表面的位置改变时，这些表达式在不损失精度的情况下进行变换。最后，在第。不及物动词我们得出结论。

二、。表面热力学及其振动膨胀

在平衡状态下，流体在给定温度下受平面无结构壁约束的自由能T型和化学势μ由提供

Ω = Ω_{b条} + γ A类,

（1）

其中Ω_b条= −光伏是巨大潜力的主要部分（P（P）作为总压和V（V）作为系统卷），A类是表面积，以及γ是表面自由能。对于这个系统，因为基底是无结构的，γ也等于流体-水表面张力。在这项工作中，表面自由能和表面张力这两个术语被用作同义词。表面吸附，即单位面积粒子数的过剩，为

Γ = {A类}^{- 1} \int_{V（V）} [ρ (第页) - ρ] d日 第页 .

（2）

在这里， $ρ (第页)$ 是密度分布和ρ表示体积流体密度。两者之间的基本关系γΓ是吉布斯吸附方程，Γ=−∂γ/∂μ.如等式。（1）系统中的平均粒子数（即V（V）),N个也分为体积贡献和表面贡献，

N个 = {N个}_{b条} + Γ A类,

(3)

具有N个_b条=ρV.数量ν是过剩体积，由关系式定义

Γ = - ρ ν .

(4)

请注意ν是多余的体积单位面积因此，有长度单位。根据等式。(4),ν如果Γ>0，则为负。与Γ一样，过剩体积与γ通过微分关系²³

ν = {(\frac{\partial γ}{\partial P（P）})}_{T型} .

(5)

对于这个系统，多余的数量γ和Γ很大程度上取决于定义体积的分割面所选的位置V（V）采用的参考区域。这同样适用于其他超额数量，例如ν然而，我们注意到P（P）和A类在这里独立于所采用的分割面约定。对于壁流HS系统，文献中经常采用两种不同的参考文献。一种是设置表面，使外部电势从无穷大降至零。其他中心围绕HS直径，σ，将表面固定在σ/2在外部电势为无穷大的区域内。此处，表面与HS颗粒表面与外墙接触的位置一致。我们分别将每个约定称为I和II，并引入符号 $\tilde{x个}$ 和 $\bar{x个}$ ⁠，用于数量x个根据I和II。因此，与这两个参考相对应的卷由

(\bar{V（V）} - \tilde{V（V）}) / A类 = σ / 2 .

(6)

中的更改γ和Γ，由等式中描述的分型面位移引起。(6)可以通过代入公式。(6)到等式。（1）和（2）分别给出

\tilde{γ} = \bar{γ} - P（P） σ / 2,

(7)

\tilde{Γ} = \bar{Γ} + ρ σ / 2 .

(8)

方程(7)和(8)混合体和表面热力学性质的转化规律。单位面积过剩体积的相应转换规律为

\tilde{ν} = \bar{ν} - σ / 2 .

(9)

等式的一个有趣特征。(9)是指 $\tilde{ν}$ 和 $\bar{ν}$ 代表纯粹的几何的移位；也就是说，它与系统的热力学状态无关。获取与γ也可以使用几何的移位，我们定义α=γ/P（P）; 然后

\tilde{α} = \bar{α} - σ / 2 .

(10)

方程(6)–(10)很容易推广。如果我们考虑分开一段距离的分割面的另外两个位置d日，关系式Eqs。(6)–(10)如果我们更换，仍然有效σ/2和d日。的值d日产生偏移，但在其他方面与ν和α从这个意义上说，我们这样说ν和α不依赖于所采用的参考，直到任意常数。使用的优势ν和α代替Γ和γ将在第节中突出显示。V B型.

对于壁流HS系统，可以方便地表达γ、Γ和其他无量纲形式的热力学量：ρ^*=ρσ^三, Γ^*= Γσ²,ν^*=ν/σ,α^*=α/σ, Ω^*=βΩ,P（P）^*=βPσ^三,γ^*=βγσ²、和μ^*=βμ，其中β= 1/k个_B类T型是逆温度（其中k个_B类是玻尔兹曼常数）。在下文中，我们去掉恒星，并假设所有的量都是以无量纲形式表示的。

受限流体的热力学性质可以用逸度的幂级数表示z（z）= Λ⁻³e（电子）^−μ，其中∧是德布罗意热长度，通过³¹

Ω = - \sum_{我 \geq 1} \frac{τ_{我}}{我!} {z（z）}^{我},

(11)

具有τ_我作为订单我簇积分。在批量中，我们有 $\frac{τ_{我}}{我!} = {b条}_{我} V（V）$ ⁠,P（P）=∑_我≥1b条_我z（z）^我、和ρ=∑_我≥1伊布_我z（z）^我.给，b条_我是的体积分量τ_我.系数b条_我是标准（体）Mayer和Montroll团簇级数理论中的团簇积分。^32–34在当前情况下，

τ_{我} = 我! V（V） {b条}_{我} - 我! A类 一_{我},

(12)

具有一_我作为的平面组件τ_我因此，除了上述体积特性方程外，我们还有

γ = \sum_{我 \geq 1} 一_{我} {z（z）}^{我}

(13)

和

Γ = \sum_{我 \geq 1} 我 一_{我} {z（z）}^{我} .

(14)

系数b条_我可以根据已知的压力维里级数系数轻松准确地进行评估。就填料分数而言η=ρ $v（v）$ ₀（带有 $v（v）$ ₀=π/6），本系列内容如下

P（P） = \sum_{我 \geq 1} {B类}_{我} η^{我} .

(15)

两者都有b条_我和B类_我表中给出了我最多八个订单。

表一。

体HS系统的簇积分和维里系数高达八阶。b条₄和B类₄都知道得一清二楚。¹¹

我	${b条}_{我} {v（v）}_{o个}^{我 - 1}$	B类_我/ $v（v）$ _o个	我	${b条}_{我} {v（v）}_{o个}^{我 - 1}$	B类_我/ $v（v）$ _o个
1	1	1	5	2168.001 435(4)	28.224 377(15)
2	−4	4	6	−22 362.204 6（2）	39.815 23(10)
3	27	10	7	243 676.514(4)	53.342(5)
4	−227.4549…	18.364 76…	8	−2 764 808.3(1)	68.528(3)

我	${b条}_{我} {v（v）}_{o个}^{我 - 1}$	B类_我/ $v（v）$ _o个	我	${b条}_{我} {v（v）}_{o个}^{我 - 1}$	B类_我/ $v（v）$ _o个
1	1	1	5	2168.001 435(4)	28.224 377(15)
2	−4	4	6	−22 362.204 6(2)	39.815 23(10)
3	27	10	7	243 676.514(4)	53.342(5)
4	−227.4549…	18.364 76…	8个	−2 764 808.3(1)	68.528(3)

另一方面，表面张力维里级数是

γ = \sum_{我 \geq 1} {W公司}_{我} η^{我} .

（16）

系数W公司_我与簇积分有关一_我通过标准反演和幂级数合成。(13)和ρ系列，低于等式。(11)]. 此外，其他表面过剩性质的级数表示，如Γ和ν在权力范围内η也可以通过通常的幂级数求逆和合成得到。因此，它们的系数很容易计算一次一_我和W公司_我已知。

三、簇积分中的曲面项

团簇积分在液体统计力学中起着关键作用。它们提供了一条准确的路径，以建立原子间势与体积流体和均匀流体的结构和热力学性质之间的关系。此外，簇积分在具有固定粒子数的系统中自然出现，并且被限制在小空腔中，即使对于小空腔N个远离热力学极限。^35,36迄今为止，团簇积分的表面积分量和表面张力维里系数是解析已知的n个= 1, 2, 3. 这里，我们精确地计算了一₄和一₅并用YSEK发布的数据进行补充，以获得比文献中更准确的结果一₆和一₇此外，这些结果使我们能够获得维里系数W公司_n个高达七阶。为了计算簇积分，我们修改了YSEK使用的方法¹⁸用于表面张力的活性膨胀。

我们的方法基于简单连接簇分解，非常方便，因为与简单连接簇兼容的配置比那些多重连接的配置更容易由MC生成，而且因为对于小型n个事实上n个=2、3、4、5、6、7，其数字分别为1、1、2、3、6、11。³⁷我们简要介绍了计算τ_n个.我们考虑n个粒子τ_n个定义为所有连接图的总和χ_n个,k个形成于n个由联接的节点（f）债券，其中（f）是迈耶函数 $（f） ({第页}_{我 j个}) = 经验 (- β ϕ_{我 j个}) - 1$ 对于粒子我和j个通过电势相互作用ϕ_ij公司每个节点表示球面坐标的积分第页_我在空间上，并由函数加权 $克 ({第页}_{我}) = 克_{我} = 经验 (- β φ_{我})$ ⁠，使用φ_我是施加在粒子上的外部电势我。则簇积分为

τ_{n个} = \int \sum_{k个} χ_{n个, k个} \prod_{我 = 1}^{n个} 克_{我} d日 {第页}^{n个} .

(17)

和往常一样，等式。(17)可以重新安排以提供较少的术语。因此，我们将所有等价的图连接起来，用一个具有多重因子的单个图来替换它们。此外，我们还确定了所有简单连接的集群，即只有一条路径连接每对节点 $χ_{n个, k个}^{(秒)}$ ⁠这些被称为树图。此外，非简单连接的簇可以看作是一个简单连接的，带有一些额外的（f）债券。因此，我们通过将其单连通部分分解为 $χ_{n个, k个} = χ_{n个, k个}^{(秒)} {F类}_{n个, k个}$ （额外（f）债券在F类_n个,k个). 我们把所有的贡献与同一个简单连接的部分结合起来，来组织被积函数，它产生

τ_{n个} = \int \sum_{k个}^{(秒)} χ_{n个, k个}^{(秒)} {F类}_{n个, k个}^{(秒)} \prod_{我 = 1}^{n个} 克_{我} d日 {第页}^{n个},

(18)

哪里 ${F类}_{n个, k个}^{(秒)}$ 是以下各项的总和F类_n个,k个相应的捐款 $χ_{n个, k个}^{(秒)}$ ⁠现在可以选择其中一个节点来扮演特殊角色并进行拆分 $\prod_{1}^{n个} 克_{我} = 克_{1} - 克_{1} {G公司}_{1}$ 具有 ${G公司}_{1} = 1 - \prod_{2}^{n个} 克_{我}$ ⁠因此，我们获得

τ_{n个} = \int 克_{1} \sum_{k个}^{(秒)} χ_{n个, k个}^{(秒)} {F类}_{n个, k个}^{(秒)} d日 {第页}^{n个} - \int 克_{1} {G公司}_{1} \sum_{k个}^{(秒)} χ_{n个, k个}^{(秒)} {F类}_{n个, k个}^{(秒)} d日 {第页}^{n个} .

(19)

等式右侧的第一项。(19)与批量贡献有关。这里一个粒子的权重是克其余的都是未加权的。第二项涉及超额非均匀贡献。这是我们在不引入分割面的情况下所能做到的。

在图的第一行。1，我们为n个=4，包含六个拓扑不同的图∑_k个χ_n个,k个在等式。(17)这里明确显示了每个图的多重因子。只有左边的两个图是树图；因此，公式。(18)仅涉及以下条款： $χ_{4,1}^{(秒)}$ 和 $χ_{4,2}^{(秒)}$ ⁠.术语 $χ_{4,1}^{(秒)} {F类}_{4,1}^{(秒)}$ 在图的第二行中绘制。1.在那里， $χ_{4,1}^{(秒)}$ 位于时代符号的左侧，并且 ${F类}_{4,1}^{(秒)}$ 在右边。我们确定了 $χ_{4,1}^{(秒)}$ 将粒子标记为1，以获得等式中的形式。(19).

图1。

查看大型下载幻灯片

方程中簇被积函数的示意图。(17)–(19).

方程中的第二项。(19)对于平面硬墙的情况可以简化。我们固定法向量 $\hat{z（z）}$ 指向流体填充区域。的功能z（z）取以下值：如果z（z）≥0，则φ(z（z）)=0和克_我(z（z）)=1，但如果z（z）<0，则φ(z（z）) = +∞和 $克_{我} (z（z）) = 0$ ⁠另一方面，G公司₁=1，如果z（z）_我<0表示至少一个我= 2, …,n个颗粒，但G公司₁否则=0。我们在z（z）=0以获得

n个! A类 一_{n个} = n个! A类 \sum_{k个}^{(秒)} 一_{n个, k个}^{(秒)} = \sum_{k个}^{(秒)} \int 克_{1} {G公司}_{1} χ_{n个, k个}^{(秒)} {F类}_{n个, k个}^{(秒)} d日 {第页}^{n个},

(20)

我们在哪里定义 $一_{n个, k个}^{(秒)}$ 作为的一部分一_n个与…有关 $χ_{n个, k个}^{(秒)}$ ⁠.

我们使用迈耶采样算法计算每个 $一_{n个, k个}^{(秒)}$ 并使用体单连通簇积分 $n个! {b条}_{n个}^{(秒)} = \int χ_{n个, k个}^{(秒)} d日 {第页}^{n个 - 1} = n个! {(2 {b条}_{2})}^{n个 - 1}$ 作为参考。此方法的第一步是固定粒子1在原点的位置。接下来n个−1个位置u个_我由Metropolis Monte Carlo生成，每个都根据分布 $π = |（f）|$ ⁠。对于给定的简单连接簇 $一_{n个, k个}^{(秒)}$ ⁠，我们使用生成的u个_我构建n个−1个粒子位置第页_我根据簇的几何形状。然后，我们评估距离 ${第页}_{我} = |{第页}_{我}|$ 和 ${第页}_{我 j个} = |{第页}_{我} - {第页}_{j个}|$ ⁠，群集函数 $χ_{n个, k个}^{(秒)} {F类}_{n个, k个}^{(秒)}$ ⁠，允许的范围z（z）粒子1的位置和长度 $我_{n个, k个}^{(秒)}$ ⁠.系数A类在等式中。(20)来自于中的集成dx公司₁第y天₁我们没有进行评估。与聚焦简单连接集群相关的贡献是

一_{n个, k个}^{(秒)} = {b条}_{n个}^{(秒)} \frac{{⟨我_{n个, k个}^{(秒)} χ_{n个, k个}^{(秒)} {F类}_{n个, k个}^{(秒)} / π^{n个 - 1}⟩}_{π^{n个 - 1}}}{{⟨χ_{n个, k个}^{(秒)} / π^{n个 - 1}⟩}_{π^{n个 - 1}}} .

(21)

使用相同的生成位置集u个_我，我们迭代地构建粒子位置第页_我对于每个 $一_{n个, k个}^{(秒)}$ 计算它对积分的贡献。我们注意到每个分母都是相同的，从而简化了计算k个，因此我们删除索引以获得

一_{n个} = \frac{{b条}_{n个}^{(秒)}}{{⟨χ_{n个}^{(秒)} / π^{n个 - 1}⟩}_{π^{n个 - 1}}} \sum_{k个}^{(秒)} {⟨我_{n个, k个}^{(秒)} χ_{n个, k个}^{(秒)} {F类}_{n个, k个}^{(秒)} / π^{n个 - 1}⟩}_{π^{n个 - 1}} .

(22)

特别是对于HS颗粒， ${b条}_{n个}^{(秒)} = {(4 π / 三)}^{n个 - 1}$ 和 ${⟨χ_{n个}^{(秒)} / π^{n个 - 1}⟩}_{π^{n个 - 1}} = {(- 1)}^{n个 - 1}$ ⁠.

四、模拟和分析

采用Rapaport算法对紧邻平面壁的硬球流体进行了分子动力学模拟。³⁸坚硬的球体被放置在相隔一段距离的两个坚硬的墙壁之间L（左）_z（z）沿着z（z）轴，具有周期性边界条件x个和年指示。在长时间的模拟运行中，我们确定了平衡流体的密度分布， $ρ (z（z）)$ ⁠, 0 ≤z（z）≤L（左）_z（z）/2，作为与墙距离的函数z（z），通过测量宽度为Δ的箱子中球体的平均数量z（z）= 0.02σ然后对两堵墙之间的结果求平均值。在50到100之间，收集每个系统的统计独立测量值，并根据这些测量值的方差确定统计误差。

密度分布在壁附近振荡，振荡衰减，密度分布值收敛到体积液体密度ρ离墙足够远。确定多余体积ν从密度剖面 $ρ (z（z）)$ ⁠，我们计算

ν (z（z）) = z（z） - \frac{{L（左）}_{z（z）} / 2 - z（z）}{\int_{z（z）}^{{L（左）}_{z（z）} / 2} ρ ({z（z）}^{'}) d日 {z（z）}^{'}} \int_{0}^{z（z）} ρ ({z（z）}^{'}) d日 {z（z）}^{'}

(23)

用于增加值z（z）直到 $ν (z（z）)$ 停止振荡，并达到相对于统计误差范围的平台。的价值 $ν (z（z）)$ 在这个平台上，被用作过剩体积的估计值。在同一时间z（z），表达式

\frac{\int_{z（z）}^{{L（左）}_{z（z）} / 2} ρ ({z（z）}^{'}) d日 {z（z）}^{'}}{{L（左）}_{z（z）} / 2 - z（z）}

(24)

用作液体体积密度的估算值ρ.

为了估计和说明由于有限的系统尺寸导致的结果中可能的系统误差，我们用L（左）_x个=L（左）_年≈ 10σ, 20σ、和40σ,L（左）_z（z）≈ 27σ和28σ对于η<3.2，以及L（左）_z（z）≈46σ和64σ对于η> 3.2. 我们没有观察到统计结果与L（左）_z（z）具体而言，关于 ${L（左）}_{z（z）}^{- 1}$ 没有显示结果对该变量的统计显著依赖性。换句话说，从具有相同L（左）_x个=L（左）_年和不同的L（左）_z（z）在所有密度的统计置信区间内是无法区分的。为了确保墙壁不会相互影响，我们使用了L（左）_z（z）远大于观测到的密度剖面振荡范围的值（与统计误差相比，振荡变得不明显的距离约为9σ在填料部分η≈0.32，约18σ最大模拟填充分数η≈ 0.49). 事实证明，这种考虑足以确保墙壁之间的任何相互作用都低于统计公差，并且L（左）_z（z）与统计误差相比微不足道，即使是在高密度情况下。

另一方面，我们确实观察到了结果对系统大小的依赖性x个−年方向，但仅在低密度下。具体来说，线性回归应用于相同结果L（左）_z（z）和三个不同的值L（左）_x个=L（左）_年揭示了线性回归模型 ${L（左）}_{x个}^{- 2}$ 因为自变量比自变量模型更适合数据 ${L（左）}_{x个}^{- 1}$ ⁠.回归线相对于 ${L（左）}_{x个}^{- 2}$ 随着密度的增加而减少，因此模型在统计上没有显示出对L（左）_x个对于η≳ 0.1. 换句话说，不同系统的结果L（左）_x个=L（左）_年在统计置信区间内无法区分η≳ 0.1.

根据这些观察结果ν这项工作中报告的结果是通过对不同的L（左）_z（z），应用线性回归模型作为函数 ${L（左）}_{x个}^{- 2}$ 并推断出 ${L（左）}_{x个}^{- 2} = 0$ 对于η<0.1，并用平均值计算结果L（左）_x个=L（左）_年≈ 20σ和40σ对于η>0.1（我们没有用L（左）_x个=L（左）_年≈ 10σ高密度）。

一次ν已计算，表面自由能通过积分方程。(5)由于在MD模拟中不直接控制压力，因此更方便对填料分数进行积分η使用修正的积分式。(5),

γ (η) = \int_{0}^{η} ν (η^{'}) (\frac{\partial P（P）}{\partial η^{'}}) d日 η^{'} .

(25)

压力对η在等式中。(25)使用Kolafa、Labík和Malijevskí（KLM）低硬球流体状态方程进行计算。^三该状态方程在整个硬球流体密度范围内精确到小数点后五位，使用该方程产生的任何误差都小于当前模拟的统计精度。在模拟中ν以等间距的一组η值（Δη= 0.005 240). 由于网格是均匀分布的，我们使用扩展的辛普森规则对等式进行数值积分。(25)。的值γ在偶数网格点，直接在所有点上使用辛普森规则计算。对于奇数网格点，使用梯形规则计算η=0.005 40，辛普森规则适用于所有后续奇数网格点。为了进一步减小数值积分误差，我们从等式中的被积函数中减去。(25)通过对的近似解析表达式进行微分而获得的被积函数γ来自白熊II（WBII）DFT，³⁹

γ = \frac{η (2 + 三 η - 2 η^{2})}{π {(1 - η)}^{2}} - \frac{自然对数 (1 - η)}{π} .

（26）

集成后γ根据方程式。（26）再加上，得出的计算值为γ通过从被积函数中减去这个解析近似，我们降低了被积函数的曲率，并提高了积分格式的数值精度。在之前的工作中，²⁵我们使用近似表达式γ但为了确保数值积分误差小于当前高分辨率数据的统计误差，需要更准确的WBII表达式。

五、结果

A.簇积分和维里系数

我们计算非均匀簇积分ã₄和ã₅与MC合作。此外，ã₆和ã₇通过补充YSEK的部分结果进行计算。¹⁸使用的精确值ã_我用于订单我=2，3，获得的数值ã_我高达我=7和b条_我（见表我)，我们计算维里系数 ${\tilde{W公司}}_{我}$ 和 ${\bar{W公司}}_{我}$ 对于我= 4, 5, 6, 7. 在表中二显示了我们使用惯例I对表面张力系列的簇积分和维里系数的结果。我们还包括YSEK的数据¹⁸用于比较。我们可以观察到本研究（PW）的结果与先前文献中的结果之间存在微小但显著的差异。请注意，其中一些小差异对应于与精确结果的比较n个= 2, 3. 值得一提的是 ${\tilde{W公司}}_{n个}$ 当我们考虑 ${\tilde{W公司}}_{2}$ 和 ${\tilde{W公司}}_{三}$ 并使用参考文献中的MC部分结果。18重新计算高阶 ${\tilde{W公司}}_{n个}$ ⁠事实上，这是第一个非实际已知的维里系数 ${\tilde{W公司}}_{4} = 13.70 (1)$ ⁠四十年前，Stecki和Sokołowski首次计算了这个积分；²¹他们获得了 ${\tilde{W公司}}_{4} = 14.045 (72)$ ⁠.

表二。

幂级数系数γ高达我=7，分格面位置采用惯例I。左边的两列对应于级数的幂z（z），而右边的对应于填充分数的幂η（维里级数）。此外， $ã_{1} = {\tilde{W公司}}_{1} = 0$ ⁠.

n个	ã_n个PW公司	ã_n个间歇泉	${\tilde{W公司}}_{n个}$ PW公司	${\tilde{W公司}}_{n个}$ YSEK公司¹⁸
2^一	−0.392 70…	0.392 69(1)	−1.432 39…	−1.432 29(7)
三^b条	2.414 53…	−2.414 6(1)	−6.097 91…	−6.100 4(14)
4	−14.360 788(6)	14.354(2)	−13.694 79(8)	−13.720(12)
5	87.609 87(3)	−87.60(3)	−22.290(2)	−22.412(8)
6	−548.818(4)	548.9(3)	−32.01(19)	−32.51(4)
7	3529.0(2)	−3530（3）	−63(17)	−65(17)

n个	ã_n个PW公司	ã_n个YSEK公司	${\tilde{W公司}}_{n个}$ PW公司	${\tilde{W公司}}_{n个}$ YSEK公司¹⁸
2^一	−0.392 70…	0.392 69(1)	−1.432 39	−1.432 29(7)
三^b条	2.414 53…	−2.414 6(1)	−6.097 91…	−6.100 4(14)
4	−14.360 788(6)	14.354（2）	−13.694 79(8)	−13.720(12)
5	87.609 87(3)	−87.60(3)	−22.290(2)	−22.412(8)
6	−548.818(4)	548.9(3)	−32.01(19)	−32.51(4)
7	3529.0(2)	−3530(3)	−63(17)	−65(17)

^一

$ã_{2} = - π / 8$ ⁠, ${\tilde{W公司}}_{2} = - 9 / 2 π$ ⁠.

^b条

$ã_{三} = 137 π / 560$ ⁠, ${\tilde{W公司}}_{三} = - 1341 / 70 π$ ⁠.

使用惯例II的表面张力的维里级数系数如表所示三我们观察到，当前工作结果的精确度提高了 ${\bar{W公司}}_{4}$ 和 ${\bar{W公司}}_{5}$ 与之前公布的数值相比，差异较大。这些微小差异可能与表中讨论的YSEK结果的微小误差有关二.精度在六阶和七阶下降。使用约定I和约定II对近似维里级数系数进行的比较表明，绝对不确定性是相同的，因为维里系数引入了误差B类_我要小得多。

表III。

表面张力系列的系数，以填充分数的幂表示。这里，使用了用于分割表面的位置的约定II。PW结果位于第二列，而第三列对应于YSEK和DLR。^18,25

n个	${\bar{W公司}}_{n个}$ PW公司	${\bar{W公司}}_{n个}$ YSEK-DLR公司
1^一	0.954 93…	0.954 93…
2	2.387 32…	2.387 32
3	3.451 39…	3.451 39…
4	3.842 27(8)	3.817(12)
5	4.662 2(3)	4.54(8)
6	6.01（19）	5.5(4)
7	−12(17)	−14(17)

n个	${\bar{W公司}}_{n个}$ PW公司	${\bar{W公司}}_{n个}$ YSEK-DLR公司
1^一	0.954 93…	0.954 93…
2	2.387 32…	2.387 32…
3	3.451 39…	3.451 39…
4	3.842 27(8)	3.817(12)
5	4.662 2(3)	4.54(8)
6	6.01(19)	5.5(4)
7	−12(17)	−14(17)

^一

${\bar{W公司}}_{1} = 三 / π$ ⁠, ${\bar{W公司}}_{2} = 15 / 2 π$ ⁠、和 ${\bar{W公司}}_{三} = 759 / 70 π$ ⁠.

B.MD和表面热力学性质的解析表达式

这里，我们展示了表面特性的模拟结果。使用MD获得的数值数据可以在补充材料此外，基于获得的维里级数系数和模拟结果，我们寻找准确描述的解析参数化表达式γ(η)和Γ(η). 在继续之前，我们应该为V（V）（因此也适用于表面特性）。自然地，一旦我们采用了一种约定并找到了适合系统表面性质的解析表达式，热力学描述就变得不那么通用了，因为考虑到描述的近似性质，当我们转换为另一种约定时，会发生一定的准确性损失。这对于γ(η)增加了一个与P（P），一个我们只知道近似值的量。这意味着，在某种约定下，一个良好的拟合函数可能会在另一种约定下转化为一个较差的拟合函数。

参考文献。31，这个问题是关于采用公约I和公约II的拟合函数的研究。在这里，为了克服这个问题，我们将重点从γ(η)和Γ(η)到数量α(η)和ν(η)通过改变原点的纵坐标，但固定其形状，可以进行简单的变换。换句话说，两者的形状α(η)和ν(η)在分割面的平移下是不变的，但对于γ(η)和Γ(η). 为了获得α（以及在不同约定之间转换），我们使用KLM-low EOS。^三

对于过剩体积，我们基于修改的维里级数提出了以下表达式：

\begin{matrix} ν & = & \frac{1}{2} - \frac{三}{2} η + \frac{339}{140} η^{2} + 一_{1} η^{三} + 一_{2} η^{4} + 一_{三} η^{5} \\ + 一_{4} η^{6} {(\frac{1}{2} - η)}^{- 1} + 一_{5} H（H） . \end{matrix}

(27)

在这里，H（H）=小时−S公司_三(小时)和 $小时 = (\frac{1}{2} - η) 自然对数 (1 - 2 η)$ ⁠、以及S公司_我(小时)是订单吗我麦克劳林系列小时[即。， ${S公司}_{三} (小时) = - η + η^{2} + \frac{2}{三} η^{三}$ ]. 方程式(27)有合适的条款η²，拟合系数有值一₁= −843/280,一₂= 111/35,一_三= −121/84,一₄=−3/1120，以及一₅= −283/1400. 中幂级数的系数η如表所示四、其中零横坐标值被忽略，因为它不相关。多余的体积ν如图所示。2其中，我们绘制了MD数据和解析表达式Eq。(27)注意，在图表的比例上，曲线与数据非常吻合。验证公式的准确性。(27)，我们在图中绘制。三两者之间的差异ν来自不同来源和当前工作MD结果。为了进行比较，我们提供了一些较旧的MD数据和现有的分析表达式（摘自参考文献。25和31。另请参见参考。40). 我们观察到MD导致PW（统计误差小于10⁻⁴)比出版的更准确，并且没有遵循现有的公式ν同时，公式。(27)在统计误差范围内很好地再现了MD结果。通过构造，拟合MD数据ν保证公式的准确性。(27)次ρ在拟合Γ时。针对以下情况详细分析了这一点α和γ.

表IV。

的维里系数ν和α（即以幂表示的级数η).

订单	ν	方程式(27)	α	方程式(28)
1	−1.5	−1.5	−0.75	−0.75
2	−2.4214…	−2.4214…	−0.192 8…	−0.1928…
3	−3.0537(2)	−3.0107	1.100 85(4)	1.1137
4	3.620(7)	3.0366	−0.372(1)	−0.52
5	−4.1(6)	−1.6022	1.6(1)	−1.0

订单	ν	方程式(27)	α	方程式(28)
1	−1.5	−1.5	−0.75	−0.75
2	−2.4214…	−2.4214…	−0.192 8…	−0.1928
3	−3.0537(2)	−3.0107	1.100 85(4)	1.1137
4	3.620(7)	3.0366	−0.372（1）	−0.52
5	−4.1(6)	−1.6022	1.6(1)	−1.0

图2。

查看大型下载幻灯片

多余体积（底部）和α比率（顶部）结果。圆圈表示MD数据，而连续线表示ν和α，等式。(27)和(28)分别是。MD数据中的误差比符号小得多，如图所示。三和4.

图3。

查看大型下载幻灯片

先前模拟结果的偏差和参数化νPW高精度MD结果。绿色的误差条表示PW中MD结果的统计误差，方块和菱形表示发布的MD数据。圆圈对应于等式。(27)、和行对应于前面的分析表达式（请参阅正文）。

对于α，我们发现基于修正维里展开的表达式

α = \frac{1}{2} - \frac{三}{4} η - \frac{27}{140} η^{2} + 一_{1} η^{三} + 一_{2} η^{4} + 一_{三} η^{5} (\frac{1}{2} - η) + H（H）,

(28)

准确描述MD数据。在这里 $H（H） = 小时 - {S公司}_{2} (小时)$ ⁠, $小时 = [一_{4} (\frac{1}{2} - η) + 一_{5} {(\frac{1}{2} - η)}^{2}] 自然对数 (1 - 2 η)$ 和 ${S公司}_{2} (小时) = - (一_{4} + 一_{5} / 2) η + (一_{4} + 三一_{5} / 2) η^{2}$ ⁠再次，我们根据订单复制条款η²，如表所示四、，而订单η^三系数的误差为1.5%。拟合参数为一₁= 227/140,一₂= −369/1400,一_三= −177/112,一₄=−11/700，以及一₅= 731/980. 在图中。2，依赖α在η已打印。MD结果与解析表达式之间的一致性(28)在绘图比例中可见。为了验证结果的质量，我们绘制了图。4等式之间的差异。(28)和MD数据，其中零对应于绘制的参考（MD结果）及其统计误差。我们可以注意到两个特性。首先，具有误差条的MD数据的高精度⪅10⁻⁵其次，从拟合表达式Eq导出的点。(28)充分描述了MD结果（在统计不确定性范围内）。在小比例尺绘制的插图中，我们比较了不同的解析表达式。这些是亨德森和普利斯克的恰当表达²⁹（惠普，虚线），乌鲁西亚³¹（U，虚线），DLR²⁵（虚线），以及Hansen Goos和Roth基于DFT-WBII的表达³⁹（HGR，点-点-灰）。此外，通过Deb直接最小化DFT-WBII得到的结果等。²⁷（DWWVB；来自表我，第3列，其中）用三角形绘制。我们注意到α来自不同来源的数据显示，偏差远大于高精度MD结果的统计误差，尤其是η>0.3之间。通过构造，表面张力（吸附）的解析表达式必须与α(ν)，即，如果Δ的表达式α产生的结果在统计误差范围内，同样适用于表面张力。

图4。

查看大型下载幻灯片

与模拟数据的偏差α误差条显示MD结果的统计误差，圆圈对应等式。(28)插图中的其他行在文本中进行了解释。

在图中。5和6显示了的结果γ=αP分别位于参考系I和参考系II中。绘制的MD数据的误差小于10⁻⁴这将使它们在绘图比例中不可见，因此不会绘制误差线。在这两种情况下，我们都观察到基于等式的公式。(28)和(10)与仿真结果吻合良好。

图5：。

查看大型下载幻灯片

MD结果和拟合函数 $\tilde{γ}$ ⁠分别以圆形和连续线表示。

图6。

查看大型下载幻灯片

MD结果和拟合函数 $\bar{γ}$ ⁠.

为了评价新的高分辨率MD结果的精度和获得的表达式的优良性γ，我们在图中绘制。7几种来源的差异γ，包括我们的解析表达式。请注意 $Δ γ = Δ \tilde{γ} = Δ \bar{γ}$ ⁠为了进行比较，DLR先前发布的MD数据^25、26也绘制了。⁴¹闭合圆基于等式。(28)。来自不同分析表达式的数据γ用线条表示。它们对应于HP，²⁹虚线，U，³¹虚线，DLR，²⁵点划线，以及HGR获得的DFT-WBII表达，³⁹分别为点-点-灰。此外，DWWVB获得了一组不同的基于DFT的结果²⁷（见表我，第3列，其中）用开圆绘制。我们注意到不同来源之间存在很大的差异，特别是η≳ 0.4. 在插图中γ在较小的范围内也很明显。

图7。

查看大型下载幻灯片

与模拟数据的偏差γ旧的MD结果取自参考文献。25和26（分别为正方形和菱形）。闭合圆由等式导出。(28)文本中解释了开放的圆和线。

我们注意到，所给出的MD结果比以前的研究精确得多，PW和参考文献的结果之间存在显著差异。25和26这种差异主要归因于高精度和小分离ΔηMD结果的ν并使用WBII表达式来减少积分误差。我们验证了公式。(28)很好地再现了MD结果。其他分析表达式不遵循新的高精度MD结果γ.

六、总结性意见

在这项工作中，研究了平面硬壁上的硬球（HS）流体系统，并给出了新的高精度结果。结果包括团簇积分[由蒙特卡罗（MC）求解]、维里级数系数和系统表面热力学性质的分子动力学（MD）模拟结果。一方面，簇积分和维里级数系数允许对系统的低密度行为进行分析描述。另一方面，为了获得HS流体整个密度范围内表面热力学的准确数值结果，模拟结果非常重要。利用维里级数系数和模拟结果，可以获得描述流体-壁界面特性的非常准确的参数表达式。我们通过参数化采用独立于参考框架的视角ν(η)（多余的表面体积）和α(η) =γ/P（P）代替吸附和表面张力（或表面自由能）。这种独立性意味着参考的改变只会在原点产生轻微的坐标偏移，从而确保拟合函数的高精度任何可以采用的参考框架。这个目标延伸到吸附和表面张力的导出表达式。

关于簇积分和维里系数，关于获得解析值的可能性的问题 ${\tilde{W公司}}_{4}$ （类似地 ${\bar{W公司}}_{4}$ 和ã₄)鉴于体积维里系数B类₄是完全已知的。精确计算 ${\tilde{W公司}}_{4}$ 这是一项艰巨的任务，因为它涉及九重积分。相比之下，我们可以提一下B类₄包含六重积分，而B类₅它也是一个九重积分，这在分析上是未知的。另一个悬而未决的问题是使用六阶和七阶系数的MC直接计算的可能性，这肯定会提高其准确性。为了将我们的方法扩展到高阶ã₆和ã₇对于简单连接的集群，必须执行。

这项对与平面硬壁接触的HS流体的研究，涉及对流体/壁热力学性质的非常准确的描述，是关于与弯曲壁接触的HS流体的未来发展的必要的第一步，特别是为了寻找其自由能中曲率项的准确描述。

补充材料

请参见补充材料用于带有MD结果的表格。

致谢

Ignacio Urrutia和Iván E.Paganini感谢阿根廷第CONICET PIP-112-2015-01-00417号拨款的支持。Brian B.Laird感谢美国国家科学基金会（批准号：CHE-1465226）的支持。Ruslan L.Davidchack在莱斯特大学授予的学习假期间完成了部分工作。分子动力学模拟是在莱斯特大学的ALICE高性能计算设施上进行的。