^{1,a)}, Andreas Greiner

^{2,b)}, Francesco Rao

^{1,c)}and Sauro Succi

^{1,3,d)}

### Abstract

We develop a three-dimensional lattice Boltzmann (LB) model accounting for directional interactions between water-like molecules, based on the so-called Ben-Naim (BN) potential [A. Ben-Naim, Molecular Theory of Water and Aqueous Solutions: Part I: Understanding Water (World Scientific Publishing Company, Year: 2010); A. Ben-Naim, “Statistical mechanics of ‘waterlike’ particles in two dimensions. I. Physical model and application of the Percus-Yevick equation,” J. Chem. Phys.54, 3682 (Year: 1971)]10.1063/1.1675414. The water-like molecules are represented by rigid tetrahedra, with two donors and two acceptors at the corners and interacting with neighboring tetrahedra, sitting on the nodes of a regular lattice. The tetrahedra are free to rotate about their centers under the drive of the torque arising from the interparticle potential. The orientations of the water molecules are evolved in time via an overdamped Langevin dynamics for the torque, which is solved by means of a quaternion technique. The resulting advection-diffusion-reaction equation for the quaternion components is solved by a LB method, acting as a dynamic minimizer for the global energy of the fluid. By adding thermal fluctuations to the torque equation, the model is shown to reproduce some microscopic features of real water, such as an average number of hydrogen bonds per molecules (HBs) between 3 and 4, in a qualitative agreement with microscopic water models. Albeit slower than a standard LB solver for ordinary fluids, the present scheme opens up potentially far-reaching scenarios for multiscale applications based on a coarse-grained representation of the water solvent.

N.M. and F.R. were supported by the Excellence Initiative of the German Federal and State Governments (DFG). N.M. acknowledges the funding from IMTEK. A.G. gratefully acknowledges funding from the Volkswagen Foundation. S.S. acknowledges financial support via the External Senior Fellow program at FRIAS. Valuable discussions with S. Melchionna are kindly acknowledged.

I. INTRODUCTION

II. THE FLUID MODEL

III. THE TETRAHEDRAL WATER-LIKE MODEL

A. The interaction potential

B. Estimating the minimum energy

IV. EQUATION OF MOTION OF THE ROTATIONAL DEGREES OF FREEDOM

A. The kinetic equation for the quaternion moments

V. NUMERICAL RESULTS

A. Simulation setup

B. From lattice to physical space-time units

C. Dynamics of hydrogen bond formation

D. Visual information

VI. LATTICE BOLTZMANN AS AN ACCELERATED MOLECULAR DYNAMICS MINIMIZER

VII. CONCLUSIONS AND OUTLOOK

### Key Topics

- Hydrogen bonding
- 15.0
- Boltzmann equations
- 10.0
- Lattice Boltzmann methods
- 10.0
- Hydrological modeling
- 7.0
- Numerical modeling
- 7.0

## Figures

The unit cell of the D3Q27 lattice. The nearest neighbors, next-nearest neighbors and next-next nearest neighbors are shown in green, blue, and red, respectively. The weighting factors, *w* _{ i }, for the D3Q27 are: *w* _{0} = 8/27 (the cell-center), *w* _{ i } = 2/27 (*i* = 1 − 6), *w* _{ i } = 1/54 (*i* = 7 − 18), and *w* _{ i } = 1/216 (*i* = 19 − 26).

The unit cell of the D3Q27 lattice. The nearest neighbors, next-nearest neighbors and next-next nearest neighbors are shown in green, blue, and red, respectively. The weighting factors, *w* _{ i }, for the D3Q27 are: *w* _{0} = 8/27 (the cell-center), *w* _{ i } = 2/27 (*i* = 1 − 6), *w* _{ i } = 1/54 (*i* = 7 − 18), and *w* _{ i } = 1/216 (*i* = 19 − 26).

Tetrahedral representation of water-like molecules sitting at lattice site **x** and its neighbor **x** _{ i }. The four arms are denoted by the corresponding normals *n* _{ k }, *k* = 1, 4 (unprimed) and *n* _{ l }, *l* = 1, 4 (primed), respectively. Blue and red code for donor and acceptor arms, respectively. Here the index *i* corresponds to *i* = 23 in Fig. 1 .

Tetrahedral representation of water-like molecules sitting at lattice site **x** and its neighbor **x** _{ i }. The four arms are denoted by the corresponding normals *n* _{ k }, *k* = 1, 4 (unprimed) and *n* _{ l }, *l* = 1, 4 (primed), respectively. Blue and red code for donor and acceptor arms, respectively. Here the index *i* corresponds to *i* = 23 in Fig. 1 .

Radial shape functions for the case σ_{ R } = 0.28 and . Note that nearest-neighbor interactions are roughly an order of magnitude weaker than next and next-next nearest neighbor ones.

Radial shape functions for the case σ_{ R } = 0.28 and . Note that nearest-neighbor interactions are roughly an order of magnitude weaker than next and next-next nearest neighbor ones.

The angular potential *V* _{ ikl } for σ_{θ} = 0.28. θ_{ ik } and θ_{ il } are shown in degrees.

The angular potential *V* _{ ikl } for σ_{θ} = 0.28. θ_{ ik } and θ_{ il } are shown in degrees.

The torque acting on the arm **n** _{ k } has two components in the direction of **e** _{ϕ}(**n** _{ k }) and **e** _{θ}(**n** _{ k }). The total torque acting on the tetrahedron is then the vector sum of the torques acting on all 4 arms.

The torque acting on the arm **n** _{ k } has two components in the direction of **e** _{ϕ}(**n** _{ k }) and **e** _{θ}(**n** _{ k }). The total torque acting on the tetrahedron is then the vector sum of the torques acting on all 4 arms.

Time evolution of the total potential energy vs the LB time steps together with the number of hydrogen bonds for four different effective temperatures, at (a) β^{−1} = 10^{−4}, (b) β^{−1} = 4 × 10^{−4}, (c) β^{−1} = 5.5 × 10^{−4}, and (d) β^{−1} = 7 × 10^{−4}. Main parameters are as follows *N* = 10^{3}, σ_{ R } = 0.28, σ_{θ} = 0.28, and .

Time evolution of the total potential energy vs the LB time steps together with the number of hydrogen bonds for four different effective temperatures, at (a) β^{−1} = 10^{−4}, (b) β^{−1} = 4 × 10^{−4}, (c) β^{−1} = 5.5 × 10^{−4}, and (d) β^{−1} = 7 × 10^{−4}. Main parameters are as follows *N* = 10^{3}, σ_{ R } = 0.28, σ_{θ} = 0.28, and .

(a) Initial random distribution of tetrahedrons and the corresponding final distribution at steady state with and without thermal fluctuations in (b) and (c), respectively.

(a) Initial random distribution of tetrahedrons and the corresponding final distribution at steady state with and without thermal fluctuations in (b) and (c), respectively.

Comparison of LB-ADR (blue) and MC (green) simulated annealing as energy minimizers. Main parameters are as follows *N* = 6^{3}, , σ_{ R } = 0.28, σ_{θ} = 0.28.

Comparison of LB-ADR (blue) and MC (green) simulated annealing as energy minimizers. Main parameters are as follows *N* = 6^{3}, , σ_{ R } = 0.28, σ_{θ} = 0.28.

## Tables

Extremal values of the angular potential , corresponding to cos θ_{ ik } = 1, 0, −1 and cos θ_{ il } = − 1, 0, 1 for the case σ_{θ} = 1.

Extremal values of the angular potential , corresponding to cos θ_{ ik } = 1, 0, −1 and cos θ_{ il } = − 1, 0, 1 for the case σ_{θ} = 1.

Article metrics loading...

Full text loading...

Commenting has been disabled for this content